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Representativeness of Wind Measurements on a 


Mesoscale Grid with Station Separations of 312m to 10000m 


by 


Steven R. Hanna and Joseph C. Chang 
Sigma Research Corp., 234 Littleton Road, Suite 2E, Westford, MA 01886 


Abstract 


A field experiment was carried out in which wind speed and direction 
were measured over flat terrain at a height of 10 m using 13 identical 
instruments spaced logarithmically along two perpendicular 10 km lines. 
Station separations ranged from 312 m to 1000 m. One-minute data from 11 


sampling periods of duration 6 to 10 hours were studied. 


The statistics showed little dependence on whether the line of instruments 
was oriented along the wind or across the wind, or whether wind speeds or wind 
directions were being analyzed. The integral time scale derived from the 
variation of the single station variances with averaging time was found to 
equal several minutes. The correlation coefficients between two stations 
separated by distance Ax were found to vary exponentially with Ax, with an 
integral distance scale on the order of 1 km. At a station separation of 
10 km, the correlation coefficient equals 0.24, 0.37, and 0.47 for averaging 
times of 1, 10, and 60 minutes respectively. These correlation coefficients 
correspend to root-mean-square differences in wind speed at the two stations of 


about 1.2, 1.1, and 1.0 m/s respectively. 


Empirical equations based on dimensional analysis are suggested for 
fitting these observations. It is found that the observations are best fit if 
two independent integral time scales are used - a boundary-layer time scale 
about 300s that best applies to small averaging times or small separations and 
a mesoscale time scale of about 1800s that applies to larger averaging times 


or large separations. 


Ls Objectives 


In order to solve most atmospheric transport and dispersion problems, it 
is necessary to assume that the wind measurements at one site represent the 
wind flow at a nearby site. The separation between these two sites can range 
from 10 m to 100 km. This assumption is sometimes quite good, but situations 
often occur where wind direction observations differ by 180° between two , 
towers in the same mesoscale network. There are a few limited studies of wind 
variability over mesoscale distances (e.g. Perry et al. (1978), Lockhart and 
Irwin (1980), Hanna (1982), and Panofsky and Dutton (1984))}, but there is much 


more theoretical and observational work needed. 


This study is part of a comprehensive research program in which the 
contributions of meteorological uncertainties to errors in air quality 
modeling and source emission estimation are being investigated. A cooperative 
two week field experiment was conducted in which 13 wind instruments were set 
up along an “L-shaped" pattern with maximum station separation of 10 km. 
Variances and spatial correlations of wind speeds and directions are 
calculated, with the objective being the development of simplified 
empirical/theoretical formulas. The accuracy of these simplified formulas is 
calculated and suggestions made for using these formulas in more 
generalized mesoscale settings. Finally, the formulas were tested with 


independent data from a wind study at Dugway Proving Ground. 
2. Previous Studies of Mesoscale Wind Variability 


The literature contains two types of studies of mesoscale wind 
variability. Both employ near-surface wind observations made by two or more 
wind monitors with separations of 10 m to 100 km. In the first type of study, 
the variances for a variety of wind monitor separations are calculated and presented, 
and a very simple empirical formula is suggested (e.g., Lockhart and Irwin (1980) 
and Hanna (1982)). In the second type of study, theoretical formulas based on 
spectral analyses are applied to the problem (e.g., Pielke and Panofsky 
(1970), Perry et al. (1978) and Panofsky and Dutton (1984)). Reviews of these 


two types of studies are given below. 


* 


Empirical Studies of St. Louis Data. During the St. Louis Regional Air 
Pollution Study (RAPS), hourly~averaged wind observations were recorded for 


one year from a network of twenty-five wind monitors. The separation of the 
wind monitors ranged from 3 km to 80 km. Lockhart and Irwin (1980) and Hanna 
(1982) report calculations of the standard deviation, o, of the difference in 
wind speed, u, and wind direction, 9, for concurrent measurements between each 
pair of stations. For example, the following procedure is used to calculate 


OM 
Oa et.) 


for any pair of stations, a and b, over the entire sampling period: 


o nn (A-1) 


where the primes indicate deviations from the average speed at each station 
and N is the. total number of time periods being analyzed. As the separation 
distance between the stations decreases, this standard deviation should 
asymptotically approach the standard deviation due to instrument errors for 


two co-located wind monitors. Lockhart and Irwin (1980) suggested the 


regression formulas: 


= 15 + 5.7lnx (A-2) 


o 
(@. 6.) 
Cen as 0.47 + 0.241nx (A-3) 
ae | 
where x is in km, 9 is in degrees, and u is in m/s. The median value of 
od is about 33° and o is about 1.2 m/s, at a median separation 
(@,-6,) (u,-u,) 


distance of about 20 km. The standard deviations of the points scattered about 
equations (A-2) and (A-3) are 0.16 m/s and a respectively. 


Hanna (1982) reported an analysis by Nappo of the same set of RAPS wind 
data, emphasizing the dependence of the standard deviations of the wind 


differences upon wind speed. Nappo defined the spatial standard deviations, 


oan and ¢ Ae using the set of concurrent measurements at the 25 monitoring 
stations. For example, the following procedure is used to calculate a for 


any given hour: 


(A~4) 


where the prime indicates a deviation from the overall average speed for the 
25 stations. The magnitude of the spatial standard deviation, Tay ts expected 
to be somewhat larger than that calculated from equation (A-2), since the 
differences in the means among the stations are included in equation (A-4). [t 
was found that the differences in wind speed and direction among the stations 
increase by a factor of three or four as wind speed decreases from 5.0 m/s to 
about 1.5 m/s. At higher wind speeds (~10 m/s), asymptotic limits of about 5° 
for o, and about 0.15 for af are evident. The following empirical formulas 


8 
fit these data: 


one (5°) + (60°/u) (A-5) 


(o,/0)* = (0.15) + (0.6/u)* (A-6) 


where u ls in m/s and Fy is in degrees. 


Theoretical Analyses of Spatial Structure 


Panofsky and Dutton (1984, Chapter 9) discuss the technical issues 
related to the differences in wind velocities measured at two points. They 
look at the problem from the point of view of the spectrum of sizes of 
turbulent eddies. They make a distinction between vertical separations, Az, 
lateral (cross-wind) separations, Ay, and longitudinal (along-wind) 
separations, 4x, between the locations of the measurements. Because turbulent 
eddies tend to maintain their identities as they travel with the wind, at a 
given separation the along-wind correlations in wind velocities between the 
measurement stations should be larger than.the vertical or cross-wind 


_ correlations. . 


A-4 


These authors ilso point out that most of the correlation between wind 
velocities at two points is due to larger turbulent eddies, with dimensions 
approximately equal to or larger than the separation between the points. By 
the same argument, the wind fluctuations in small eddies at one point are 
uncorrelated with the wind fluctuations ir the same size of eddies at a 


distant point. 


Because of the fact that the larger eddies move with the wind flow and 
Slowly change with time, the maximum correlation between two time series of 
wind velocities at two points separated by an along-wind distance, Ax, may 
occur at a time lag of about A4x/U, where u is the mean wind speed. This is 
the approximate time required for an air parcel to travel from the first point 


to the second point. 


- The issues raised in the previous three paragraphs have led Panofsky and 
Dutton (1984) and others to the hypothesis that the spatial correlations are 
functions of the eddy size, the wind speed, and the separation between “he 
wind observation points. They use a mathematical expression known as the 
coherence, which “is a measure of the square of the correlation between the 
Fourlér component of two time series with their phases adjusted to obtain 
maximum correlation." The coherence is a function of eddy frequency, n, and 


is given by: 


2 2 
coh(n) = ico’(n) + Qitn)} (A-7 ) 
3, (n) S,(n) 


where S, (a) and S, (n) are estimates of the spectral density of the two time 
series at points 1 and 2, Co(n) is the cospectrum, and Q(n) the quadrature 
spectrum. The frequency, mn, can have units of cycles per second or radians 
per second. The cospectrum refers to the correlation due to in-phase 
fluctuations and the juadrature spectrum refers to the correlation due to 
fluctuations that are 90° out of phase. The coherence ranges from zero, for 
no correlation, to one for perfect correlation. In order to calculate coh(n) 
for two long time series, speciaiized computer software fer carrying out time 


series analysis is needed. 


Pielke and Panofsky (1970) generalized an empirical expression for the 


coherence, coh(n), based on a suggestion by Davenport (1961): 
coh(n) = exp(-a, nAx, /u) (A-8) 


where the subscript i refers to the direction component and the dimensionless 
“constant”, ay is called the “decay parameter." Observations show that "a," 
is typically in the range from 1 to 10. The dimensionless variable, nx, /u, 
is often referred to as the reduced frequency. It is found that the 
“constant,” ay is in fact a function of station separation, Ax. Other 


investigators have found a, to be a function of the ambient stability, the 


i 
turbulence intensity, and the integral time or distance scale of the 
turbulence. For example, Perry et al (1978) propose the following formula for 


horizontal station separations, Ax: 


a = (65¢ /u) + (6.3¢ Ax/uL ) (A-9) 
W v y 


where e. and o, are the standard deviations of turbulent velocity fluctuations 
in the vertical and cross-wind directions, and Ly is the lateral integral 


distance scale of the lateral turbulence. 
3. Proposed Formula 


The relation between wind speed, u, or wind direction, 6, observed at two 
stations at positions 1 and 2 can be expressed in terms of a variance or a 
correlation. The relationship between the two quantities can be derived by 
beginning with the identity: 


+ u - 2 u,’u.’ (A-10) 


where the nrimes refer to fluctuations (i.e., the means have been already 


subtracted) and the overbar refers to a time average. If the two variances 


se and a are equal, then equation (A-10) can be written in the form: 


u 
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Tru 20, =1- Rio (A-11) 
where oe = u’? and the correlation coefficient Rio = uy /uy’ 7a. Thus if the 
correlation equals 1, 0, or -1, the ratio Caye/2e will equal 0, 1, or 2, 


. respectively. For two stations that are separated by such a large distance 


that the correlation equals zero, then Equation (A-11) states that es = an *. 


After investigating the semi-theoretical formulas (A-8) and (A-9), it was 
decided that there was so much adjustment of “constants” taking place that it 
was best to begin with a simpler formula based on dimensional analysis. 

This simpler formula does not explicitly include the frequency, n, but does 
account for frequency effects through the implicit assumption that wind speed 
autocorrelograms are exponential (i.e., R = exp(-AtvT,), implying that a 
so-called Markov-energy spectrum is valid) and are completely determined once 


an integral time scale, T or space scale, Ay is specified. The following 


[’ 
dimensionless relationship can then be postulated: 


2 2. 


Au a) (A~12) 


T’ 


where f is a wpniversal dimensionless function, and variables and parameters 
are defined in the following way: 


a Cu, - uy ee where us and uy’ are wind speed fluctuations 
at positions 1 and 2. 
2 2 2 : oe 
Fass ie. aa 2 )/2 is the average turbulent energy at positions 
1 2 


1 and 2. 


4x is the horizontal separation between positions 1 and 2. 


Ty is the integral time scale of the time series of wind speed 


fluctuations measured at position 1 or position 2. 


Ay is the integral space scale of the correlogram of wind fluctuation 
differences (u,’- uy") measured between positions 1 and i for 


several positions. 


qT, is the averaging time. 


A-7 


The averaging time, tT. is very important because it is typically of the 
same order of magnitude as the integral time scale (i.e., within the range 
from 1 min to 1 hr). As a approaches zero, the ratio Ou 20° reaches a 
maximum at any given station separation, since there are poor correlations 


between wind fluctuations in smaller eddies at two positions. 


Exponential shapes are proposed for all correlation functions: 


R(Ax/A,) exp (-Ax/A, ) (A-13a) 


R(At/T,) = exp(-at/T,) (A-13b) 


The following approximation to Taylor’s equation is used to account for the 


effects of averaging time on variances: 


o°(T,)/e*(0) = (1 + T,/2T)) | (A-14) 


For a given separation, Ax, the spatial correlation R(Ax/A, ) is also going to 
depend on averaging time, tT. since the effects of the fluctuations due to 
small eddies will drop to zero at large averaging times. Consequently, 
R(AX/A, ) should approach unity as qT, increases, and the second term in the 


following empirical formula is proposed to account for this effect: 


R(T, Ax) =(u, uy Ce v, I(T, Ax) 


2 


= exp(-Ax/A,) + (1-exp (-Ax/A,) )g(T uA, ) (A-15) 


where the dimensionless function a(T uA, ) approaches zero as T approaches 


zero and approaches unity as T, becomes very large. 
4. Field Experiment at Hereford, Colorado 
The research programs summarized in Section 2 all took advantage of wind 


observations made at two or more locations. However, none of these wind 


observations satisfied all of the following desired characteristics of a study 


in which differences in wind observations are being calculated. 


° Flat terrain with few obstructions. 

. All wind instruments alike, with relatively fast response and 
adequate calibration procedures. 

° Five or more wind stations along a line with spacings ranging from a 
few hundred meters to 10000 meters. 

. Two perpendicular lines of wind stations. 


. Several days or more of one-minute averaged wind observations. 


A field experiment was designed to satisfy these characteristics, using. 
instruments and technicians from the National Center for Atmospheric Research 
(NCAR). Thirteen identical wind stations from NCAR’s Portable Automated 
Mesonetwork (PAM) were set up in flat, open farmland (mostly covered by short, 
winter-weathered grass) near Hereford, Colorado, and operated for the two week 
period between 30 March 1990 and 14 April 1990. The instruments were 
installed, maintained and calibrated by experienced NCAR technicians. The 
anemometers were located at an elevation of 10 m above the ground. A 
schematic diagrams of the relative instrument locations is given in Figure 
A-1. The slope of the terrain is less than about 1% over the entire network. 
It its seen that logarithmic spacings (312.5m, 625m, 1250m, 2500m, S00Om, 
10000m) are used for the instruments along each leg of the “L"-shaped network. 
The perpendicular axes were used in order to determine if there was a 


significant difference between the along-wind and cross-wind statistics. ~ 


The resulting "Hereford" dataset contained values of N-S and E-W 
components of the wind velocity at one minute resolution. First, these wind 
components were converted to wind speed, u, and wind direction, 9, for each 
minute. Because the threshold of these instruments was about 0.5m/s, periods 
with reported wind speeds less than 0.5m/s were flagged and not used in the 
analysis. S-minute averages of u and @ were made and the resulting time 
series for each station plotted for the entire 16 day duration of the study 
(see Figure A-2 for an example of the time series plots for 31 March 1990 
Station 6). Visual inspection of the 15 days of time series in the figure 
revealed that, most of the time, the wind speed and direction are quite 
variable and unsteady in time. For example, during the afternoon on 1 April 
and 10 April, the wind speed went through a cycle in which it increased from 
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Figure A-1. Schematic diagram of wind station locations at the Hereford site. 
Terrain sloped slightly downward from west to east with an 
average slope of 1%. North is towards the top of the figure and 
there is 10 km spacing between stations 1 and 7 or stations 7 and 
13. 


A~-10 


nearly zero to about 10 m/s and decreases back down to zero again. This large 
half-sine-wave would totally dominate any statistical analysis of these wind 
data. 


Because it became clear that any time series analysis of the 16 days 
would be overly influenced by diurnal changes and synoptic effects, it was 
decided to confine the analysis to time periods on the order of 6 £5 10 hours, 
when the wind speed and direction appeared to be relatively steady. The 
eleven steady-state periods that were chosen are listed in Table A-1. Two of 
these periods are indicated by thick lines in Figure A-2. Median values of 
wind speed, u, wind direction, 9, standard deviation of wind speed, oy and 


standard deviation of wind direction, ¢ are also listed on the table. These 


9’ 
eleven periods are seen to cover a wide range of wind speeds (2.5 to 10.5 


m/s), wind directions(170 to 330°) and times of the day. 
5. Analysis of Hereford Data 


This section presents the results of the analysis of the Hereford data 
and the next section presents some empirical formulas that fit this dataset. 
As a first step, variances and time and space correlations were calculated for, 
averaging times of 1, 10, and 60 minutes. A linear trend (estimated from the 
data from all 13 stations) was removed from the data from each run. Because 
steady-state periods have been selected and linear trends have been removed, 
much of the influence of larger mesoscale and regional eddies has been removed 
there statistics. As will be shown, mesoscale eddies with time scales 


less than the sampling time are still strongly reflected in the statistics. 
5.1 Single.Station Variances as a Function of Averaging Time 


The variances in the wind speed and wind direction time series were 
calculated for each of the eleven “steady-state” runs, for each of the 
thirteen stations, for averaging times of 1, 10, and 60 min. The data are 
presented in a Table A-2, in the form of ratios of the variance for 10 or 60 
minute time periods to the variance for the one minute time period. The 


medians over the thirteen stations are listed. The overall medians at the 


“Anil 


Table A-1 


Steady-State Periods from the Hereford Dataset Selected for Analysis. 


Meteorological Parameters Represent Medians over 13 Stations 


Wind Speed Wind Direction 


Date/Run Time — u(m/s) aC )} & (m/s) a, ( ) 
0330 00-10 2.5 130 0.7 17 
0330 14-24 3.0 215 0.8 17 
0331 00-08 3.5 200 0.45 10 
0331 16-24 7.5 325 1.2 13 
0404 06-12 5.5 90 1.2 13 
0404 12-20 7.5 135 1.2 10 
0406 18-24 4.0 215 0.9 14 
0409 14-24 10.5 10 1.3 8 
0411 06-16 3.0 170 0.9 13 
0412 00-08 .5 140 1.1 11 
0413 14-20 .5 330 1.1 10 
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Figure A-2. 
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Time series of five-minute averaged wind speed (u) and wind 
direction (@) for Station 3 for 31 March 1990 at the Hereford 


site. Thick lines indicate steady-state periods selected for 


analysis. 


AW13 


bottom of the table suggest that there is little difference between the 
statistics for wind speed and wind direction, with medians of the variance 
ratios of about 0.68 for o7(10 min)/e“(1 min) and about 0.39 for (60 
min)/o"(1 min). The run-to-run variability in the median variance ratios is 


about + 0.20, with no dependence on wind speed, wind direction, Ty Ty OF 


9’ 
time of day. For a given day, the station-to-station variability in the 


variance ratios is typically about + 0.05 to 0.10. 


Table A-2 suggests that observed median values of (10 min)/o7(1 min) 
and o*(60 min)/o-(1 min) are 0.68 and 0.37, respectively. Solving Equation 
(14) for the integral time scale Th we obtain T, = 9 min for qT. = 10 min 
and Ty = 15 min. for Tt, = 60min. This trend is duplicated in all types of 
calculations in this analysis, i.e., the derived time and distance integral 
scales increase as the averaging time, tT) or station separation, Ax, 
increase. Persistent mesoscale eddies are associated with turbulence at time 
scales of several hours and distance scales of several tens of kilometers. 
Superimposed on these “baseline mescscale eddies" are the smaller turbulent 
eddies normally thought to be associated with the atmospheric boundary layer. 
Consequently it is assumed that the turbulent velocity field is described by 
two times scales, one (T, 1) associated with boundary layer turbulence, and 
another (Tro) associated one eer ee eddies. These two time scales are 
assumed to contribute equally to¢e mT, ), leading to the following modification 
to Equation (A-14): 


2 
ee 0.5 ; 0.5 nena 
o7(0) ad ee a eer 


— ee dataset is best-fit by Thy # 300s and Tho # 1800s, which yield 
- S(10)7e" (1) = 0.71 (slightly above the median observation of 0. 68) and 


gv 2 (60) /e* (1) = Q.35 (slightly below the median observation of 0.37). 
5.2 Statistics for Pairs of Observing Stations 
The spatial statistics to be presented are all keyed to wind station 


(7) at the corner of the “L" of the network shown in Figure A-i. The results 


in this section are given for the spatial correlation coefficient, 
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Run 
Date/Time 


0330/00-10 
0330/14-24 
0331/00-08 


0331/16-24 
0404/06-12 . 


0404/12-20 
0406/18-24 


0409/14-24 


0411/06-16 


0412/00-08 | 
0413/14-20 | 


Overall 
Median 


Table A-2 


Ratios of Variances o?(T Ie" (1 min) for Averaging Times, qT, of 10 and 

60 min., for Wind Speed and Wind Direction at the Hereford Site. Results 

are Given for Each of Eleven Runs. The Median of the Ratios for the Thirteen 
Stations is Given. The Scatter of the Ratios for the Thirteen Stations about 
any Median has a Standard Deviation of about 0.05 to 0.10. Medians over all 


Dates are Given at the Bottom. 


2 2 P 2 2 
o, ear (1 min) an (Tle, {1 min) 
Wind Speed : Wind Direction 
10 Minute 60 Minute 10 Minute 60 Minute 
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R(4x, T_) =u “Wy Ao ,o  , which are functions of the separation, Ax, of 
a i uy uw 


stations i and 7, the averaging time, Te and the integral scale of the 
turbulence. Averaging times of 1, 10 and 60 min are used. Results are 
presented separately for the E-W and N-S legs of the “L". Rather than showing 
the correlation for all eleven time periods or runs, we present the median of 
the 11 periods. These correlations, are plotted in Figure A-3. There does 
not appear to be a strong dependence on E-W or N-S leg, or on whether wind 


speed or wind direction is being plotted. 


As expected, the calculated correlations are lowest for the shorter 
averaging times and the largest station separations, due to the fact that the 
dominant turbulent eddies are characterized by space scales of a few hundred 


meters and time scales on the order of a few minutes. 


At the 10 km separation distance, the correlations average about 0.5. 


2 


Because of the relation between spatial correlation, 2A Spatial variance, oc 


LZ Au 
and variance, °°, indicated by Equation (A-11), it follows that, at 10 km 
separation, ts = a. Knowing that oT, = 1 min) = 1.2m*/s*, and using 
the median values of oT vert min) in Table A~2, it can be concluded that 
Oyu ~ 1. tm?/s* for T, = 1 min, o,.* ~ 0.8m°/s* for T, = 10 min., and 
ey 0. 4m/s* for ay = 60 min at the 10 km separation distance. This value 
of Ri Q.7m/s for one-hour averages is close to that found with the St. 


Louis RAPS data, as discussed in Section 2. 


At small averaging times, tT, Equation (A-15) suggests that the correlation 
coefficient for wind speed fluctuations at two points separated by a distance, 
Ax, should be an exponential function of Ax. Therefore, if In R is plotted 
against Ax, a straight line should be seen. Observations of R from this wind 
hetwork are plotted in Figure A-4, for the six Ax values (312.5, 625, 1250, 
2500, S000, and 10000 m) and three averaging times (1, 10, and 60 min) of 
interest. Because the correlations in Figure A-3 did not indicate nuch 
dependence on E-W or N-S direction, and were similar for wind speed and wind 
direction, all runs, legs, and variables are combined in this table and 


figure. 
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Figure A=3. 


a) 


2000 20000 200 2000 20000 
E-W ax (m) 


Spatial correlation coefficient, R, as a function of station 


separation for E-W and N-S legs of Hereford monitoring network. 
Wind speed (u) correlations are on the left, and wind direction , 
(@) correlations are on the right. Medians over all eleven runs 


are shown for averaging times of 1 min. (Long Dashes), 10 min. 


(Short Dashes), and 60 Min. (Solid Line). The standard deviation 


of the scatter of the 11 points about each line is about 0.2. 
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Figure A-4. 
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Spatial correlation coefficients, R, for various spatial 
separations, Ax, and averaging times, Ti (a 1 min; * 10 min; + 60 
min), from Hereford site. Plotted are the medians over eleven 
runs. E-W and N-S legs, and wind speed and wind direction 
observations. Typical scatter of all the data about each line is 
about 0.2 at a correlation of 0.5. Empirical curves from 
Equation (A-17) are drawn (dotted, Ta = 1 min.; dashed, Ta = 10 
min.; solid, Ta = 60 min. ). 
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The points in the figure follow a straight line only for the largest 
averaging time (60 min). For the smallest averaging time (1 min), there is a 
much more rapid drop-orf in correlation at small separation cistances. This 
behavior can be explained by the same two-scale phenomenon discussed under 
Section 5.1 in the analysis of the effects of averaging time; i.e., the 
results are influenced by a smaller eddy scale, Ay representative of boundary 


layer processes, and a larger eddy scale, A representative of mesoscale 


2° 
fluctuations. However, the effects of the smaller scales tend to disappear at 


large averaging times, Ts when the condition T > A/a is satisfied. Thus, 


equation (A-15) appears to be satisfied, with the limits: 


~Ax/A 
e 


R( Ax, rt >» 60 nin) 2 


R(ax, T, >t min) = 0.5e 77") + 0,5e 42 
The following empirical equation fits these data and has the proper 


asymptotic behavior: 


- oA) (1 + (Tuvan, )* (A-i7) 


where Ay = 300 on, A = 12000 m, and a = S when u = 5 m/s. The predicted 
curves are drawn on Figure A-4. The value of a = 5 is consistent with 
“Pasquill’s Beta" = 4, which is the known proportionality factor between 
Lagrangian and Eulerian scales. Because the actual eddy diameters are equal 

to about five times the integral scales A, the smaller eddy diameters would de 
about 1500 m and the larger eddy diameters would be about 50 xm. Of course the 
confidence limits on the data (+ about 0.2) are so broad that a much wider 
range of scales, Ay and As, is possible. Furthermore, these results are 
limited to this particular site and time of year, and to the ranges of Ax, tT 
and meteorological corditions that we have considered. The sampling time (a 


maximum of 10 hours) also imposes a maximum limit on integral scales. 
6. Recommended Empirical Formula for Variances 
Our analysis of the Hereford data has suggested that the variances and 


correlations for wind speed and wind direction are similar. In addition, we 
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have found that, despite the theoretical prediction that along-wind 
correlations will be larger than cross-wind correlations, there is no clear 
dependency of the correlations on wind direction (or any other meteorological 
variable). This apparent lack of dependency may be due to the fact that this 
effect is overwhelmed by the natural variability in the observations. 


The similarity relations in Equation (A-8) have been shown to be valid, 
resulting in Equation (A-16) for the effects of averaging time and Equation 
(A-17) for the effects of station separation. These equations can be combined 


into the following general equation: 


2 

Faux Th) 0.5  , 05 

ett wi), [PP a, alae 

uoa 
sibs tlhe eM, 
Le +o.5 ¢ ae) (A-18) 

(1 + (Tuan, )*) 
a 1 
where Thy ~ 300 s Tr ~ 1800 s 
A, ~ 300 m Ay ~ 1200 m 


a (Lagrangian-Eulerian scale) ~ 5 


This equation should be tested with independent data from a different 
Site. Slightly different values of the time and distance scales may be 
appropriate at a different site. For any given time and place, the confidence 
limits on the results would be expected to be in the same range (about 20% to 
30%) has been found here. 


7; Test of General Equation with Independent Data 


A set of independent wind data from the so-called XM21 field study at 
Dugway Proving Ground, Utah, were provided by C. Biltoft. Data were available 
from two towers, separated by 500m, on a relatively flat test range. 
Instruments were located at heights of 2, 4, 8, 16, and 32m, and measurements 
were made at a frequency of one per second. Because the towers were 


constructed of scaffolding, it was necessary to disregard periods when the 
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wind blew through one of the towers. Two periods of valid data were 
analyzed--from 1110 to 1651 on 11 April 1989, and from 0937 to 2059 on 2 May 
1989. Both periods were marked by wind speeds of about 3.5 m/s and relatively 
strong turbulence intensities (o 7a ~ 0.3 and or as 30°). 

Because the integral scales of the horizontal components of turbulence 
are not strongly dependent on height, the results from all five levels are 
combined in our analysis. This assumption will cause a slight error, since 
the various ratios 0? (T,)/e*(1sec) and correlations are observed to increase by 
about 10% to 20% between heights of 2m and 32m. In addition, the wind speed 
and direction results appear similar and are combined in our analysis. The 
combined data results (two variables, two days, five levels) are presented in 
Table A-3. 


The predictions of o° (7, )/o* (18sec) in the table have been made using 
Equation (A-16), with the same time scales derived from the NCAR wind network 
(T, = 300s and Ts = 1800s). There is a + 10% agreement between the 
predictions and the observed medians, and the predictions are within the + 
standard deviation error bounds of the observations. It can be concluded that 
the variances o°(T,) show similar behavior at the two sites, and that Equation 


(A-16) can adequately simulate this behavior. 


Equation (A-17) (for spatial correlations) does not transfer as well to the 
new site. Its predictions are shown under column (2), assuming the parameters 
(Ay = 300m, Ay = 12000m, and a = S) derived from the NCAR wind network. 
Clearly the predictions of Ray ita) are too low by about 10 to 30%. If the Ay 
A, and "a" parameters are "tuned" with these new data, the predictions under 
column (3) are all within + 4% of the observations. However, the smaller of 


the length scales, A,, has to be increased from 300m to 1000m, and the 


Lagrangian-Eulerian eee "a" has to be reduced from 5 to 0.2. It appears 
that the weakest part of Equation (A-17) is the correction term, (1 + 

(ut 7aa,)7)7*, for averaging time, Ts: This term is intended to cause an 
increase in the spatial correlation as averaging time increases. However, the 
functional form for this correction term is not obvious and more thought is 


clearly needed. 
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Table A-3 


Observations and Model Predictions for Dugway Proving Ground Wind Data, 
with Tower Separation of 500m, for All Levels Combined, 
and Wind Speed and Direction Combined. Medians are Listed for the 
Observed Values. Model Parameters are given in the 


Footnotes. 
o7(T_) o7(T.) 
a a 
Averaging o*(1sec) o*(1sec) Ry (T_) Ry (T_) 
Time ua ua 
tT. Gbserved Predicted Observed Predicted 


(2) (3) 


0.78 0.57 0.77 


60 sec 0.83 , 0.58 0.86 


600 sec 


Footnotes: (1) Equation (A-16) is used, with the same values for the 


parameters as derived from the NCAR data (T, = 300s and 


1 
TS = 1800s) 


(2) Equation (A-17) is used, with the same values for the 
parameters as derived from the NCAR data (A, = 300m, 


Ay = 12000m, and a = 5) 


(3) Equation (A-17) is used, with A, = 1000m, Ay = 12000m, and 
a =0.2. 
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Definition of Terms 


A,B,C,D coefficients which are determined by fitting the observed 


wind tunnel experiments to a curve 


COV covariance of 
xy xy 
a 
cv2 squared coefficient of variation of x, = 
x 
2 
; c 
oy squared coefficient of variation of y, 4 
y 
h relative humidity 
M, initial mass 
My mass of zinc or phosphorus aerosolized, as measured in the 
wind tunnel experiments 
M,, 
MYF munition yield factor, the ratio x 
° 
Q, mass emission rate 
Q 
q a 
° 
x mean of x 
y mean of y 
YF yield factor, theoretical adjustment to mass of aerosol based on the 


ability of the active ingredient to absorb water vapor 


AM mass lost during burn, initial mass less final mass, as measured 
by the load cell 


variance’ of the product xy 
variance of x 


Ce variance of y 


Appendix B-1 
UNCERTAINTY ASSOCIATED WITH EMISSION RATE ESTIMATION 


dy Introduction 


To estimate the uncertainty of an entire model it is necessary to 
evaluate the uncertainty of its components. The input data for atmospheric 
dispersion models consists of emission rate data and meteorological data. In 
this discussion we address the issue of uncertainty associated with the 
estimation of emission rates of aerosols from smoke munitions. Emission 
rates are generally expressed in terms of weight per unit time; in this case 
they are expressed as grams of active ingredient, (e.g., phosphorus) per 


second. 


The issue of estimating the emission rate of obscurant munitions is 
complex, aside from the issue of estimating the associated uncertainty. There 
are several reasons for this. First, the critical property which provides the 
obscurant effect is not directly emitted by the munition; it results from an 
interaction of the active ingredient, e.g., red phosphorus, and moisture in the 
ambient atmosphere to form a dense smoke cloud. Second, munitions contain 
other materials which burn simultaneously but do not contribute directly to the 
obscurant effect; thus measuring total weight loss over burn time only 
indirectly measures the amount of active ingredient which has been released. 
Third, to experimentally determine the amount of active ingredient released, 
the mass of the active ingredient in the entire smoke cloud must be determined 
and this measurement can only be carried out in a wind tunnel. Thus, data from 
the wind tunnel experiments are extrapolated to the field setting. These 
issues affect emission rate estimation as well as contribute to the associated 


uncertainty. 


This discussion covers the following topics: how emissions are 
measured in both wind tunnel and field experiments, the model used for 
estimating emission rates, a technique for estimating uncertainty (using 
existing data as an example), and conclusions and recommendations for the 
collection of additional data which would enhance the uncertainty analysis. 
Several reports were reviewed for this analysis, but the infermation on the 


method of measuring and modeling emission rates is based on two reports: Basic 


Snoke Characterization Test (DPG-FR-77-311) (DPG,1978) and Methodology 
Investigation Final Report Validation of a Transport and Dispersion Model for 
Smoke (DPG-FR-702) (Carter et al., 1979). These reports contain data on 
three different types of smoke obscurant munitions: white phosphorus, red 
phosphorus, and zinc oxide-hexachloroethane-aluminum. Because only limited 
amounts of appropriate data are available for the uncertainty analysis, the 


conclusions drawn from this analysis must be considered tentative. 
2. Measurement of Smoke Munition Emissions 
2.1 Wind Tunnel Experiments 


The most detailed measurement of emission rates of submunitions were 
conducted in wind tunnel experiments. Although field measurements have been 
collected, these data are not adequate for developing a model. Thus, the 
model used to predict emission rates is based on the wind tunnel 
experiments. The experiments and the data described in this section are 


summarized from DPG (1978). 


The wind tunnel experiments were conducted at the Dugway wind tunnel. In 
the experiments, single submunitions (e.g., an individual component which 
contains the smoke producing compound) of zine and white phosphorus were 
burned. [In actual field use of smoke munitions, multiple submunitions are 
loaded into a canister and burn simultaneously. For red phosphorus, three 
submunitions (wedges) were burned simultaneously. The wind tunnel tests were 
conducted at a wind speed of 2 m/s and the ambient temperature, relative 
humidity, (“RH), and barometric pressures were recorded. Pre-weighed 
submunitions were placed on a load cell in the center line of the tunnel and 
fired remotely. The load cell recorded the elapsed time and weight loss of the 
submunition as it burned. The end of the tunnel had a barrier that contained a 
grid of holes through which the air in the tunnel passed. The air passed into 
sampling lines and through collection devices, impingers. (See Figure B-1.) The 
impingers were then analyzed for zine or phosphorus, the appropriate active 
ingredient. The concentrations in all impingers were summed to determine the 
total amount of active ingredient which had been aerosolized during the burn. 
Data collected on each test included: weight loss with burn time, total zine or 


phesphorus, tunnel wind speed, temperature, and relative humidity. 
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Figure B-1. Ofiagram of Wind Tunnel. Source: DPG, 1978 


Transmittance, which is routinely measured in the field, cannot be measured in 
the wind-tunnel because the concentrations are too high. The types of 


submunitions tested (as reported in DPG (1978)) were: 


155mm HC M1 Canister (zinc) 

155mm HC M2 Canister (zinc) 

10Smm HC M1 Canister (zinc) 

6 inch WP Wick (white phosphorus) 

2.75 inch Rocket WP Wick (white phosphorus) 
8imm Navy RP Wedge (red phosphorus) 

155mm Navy RP Wedge (red phosphorus) 

81mm German RP Wedge (red phosphorus). 


For each type of submunition, two burns were conducted. 
2.2 Field Experiments 


The field experiments were conducted at the Horizontal Grid, Dugway 
Proving Ground, Utah. The experiments described in this section are summarized 
from DPG (1978). The data were collected by means of photography, aerosol 
sampling with aerosol photometers, particle size analyzers, and impingers. 
Samplers were located 1.22 meters above the ground. Motion pictures recorded 
the size and shape of the smoke cloud; aerosol photometers recorded total 
particle concentration; particle size analyzers recorded size distributions; 
and impingers measured zinc or phosphorus concentrations. The mass of zine or 
phosphorus from these impingers cannot be summed to calculate total Zinc or 
phosphorus aerosolized as was done in the wind tunnel experiments because the 


total smoke cloud is not collected by the impingers. 


Figure B-2 (DPG, 1978) shows the layout of the instruments for these 
experiments. The sampling line was always oriented perpendicular (+ 45°) to 
the prevailing wind direction. The impingers were located halfway between the 
aerosol photometers. As in the wind tunnel experiments, the munitions were 
placed in a load cell so that data on weight loss with burn time were recorded. 
The types of submunitions tested were the same as those used in the wind tunnel 


experiments. For each type of submunition, two burns were conducted. 


Figure B-2. - 
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Grid for munition cloud characterization. Source: 


ae i LL, 


1978 


The meteorological data recorded in these experiments included: wind 


speed and direction, temperature, and relative humidity. 
3. The Emission Rate Expression 


The critical parameter to be derived from these experiments is the 
amount of smoke generated per second of burn time. In its simplest 
form, this would be calculated by dividing the mass lost by the burn time 
to derive a value in units of grams per second. However, the amount of smoke 
generated is more than a function of the mass burned; it is also a function 
of the interaction of the zine or phosphorus with water vapor to create the 
smoke. Zinc and phosphorus are hydroscopic, absorbing upwards of four times 
their mass in water vapor. Thus, the actual emission rate of interest is 
total amount of hydrated smoke generated per second. To add one further 
complexity to the equation, the burn rate of the munitions is not constant 
with time. There is an initial growth >eriod, followed by a plateau of 
relatively constant rate of emissions, - 2llowed by a rapid decline. In 
field tests, there are data which indicate munitions can burn unevenly; this 
is characterized as "flashing" toward the end of the burn time by the 
observers. | 


The following discussion covers two topics. The first topic presents 
- the parameters in the emission rate expression and the second discusses the 
uncertainty associated with those parameters anda presents a technique for 


computing the uncertainty in the overall emission rate value. 
3.1 Input Parameters 


The emission rate of hydrated smoke is a function of the initial mass 
of the munition, the fraction of this mass which is the aerosolized zinc or 
phosphorus, the increase in mass of the aerosol due to hydration, and burn 
time. The expression which was developed from the wind tunnel experiments to 
estimate the emission rate is as follows (DPG, 1978): 
2 | 


3 
Q. = Mo MYF YF (A/t, + 2Bt/ ty + 3Ct 7b, + 4Dt ty, 


-) (B-1) 


where: Q, = mass emission rate 
M, = initial mass 
MYF = munition yield factor 


YF = yield factor ; 
A,B,C,D = coefficients which are determined by fitting the 
observed data in the wind tunnel experiments to a curve. 


They are constrained to sum to unity. 


The munition yield factor, MYF, is the ratio of the mass of zinc or phosphorus 
burned, Me to the initial mass of the munition, Mo: MS is determined from the 
wind tunnel experiments and is calculated by summing the mass of 2ince or 
phosphorus collected by all the impingers. These wind tunnel values of M.. are 
used to predict the emission rate, Q,. for field experiments, as MY cannot be 


measured in the fleld setting. 


The yield factor takes into account the hydroscopic growth of the 
aerosolized zinc or phosphorus and is a function of ambient relative humidity. 
For the three different types of smoke munitions, yield factors are based on 
theoretical calculations or experimental data. The conceptual definition of 
the yield factor is: 


yF = mass of smoke (including aerosol) 
mass of starting material 


In the case of phosphorus, the final material is primarily hydrated 
orthophosphoric acid, (HPO, + nH,0). According to DPG (1977), it is 
assumed that the hydration of the droplets of HPO, proceeds until their 
aqueous vapor pressure equals the partial pressure of He in the atmosphere. 
The mass of hydrated HPO, is then constant (assuming no other environmental 


changes occur) and the yield factor can be expressed as: 


mass HPO, + nO 


YF = mass P 


DPG (1977) contains the following example of how a yield factor is 
calculated. Yield factors for phosphorus can be readily computed from tables 


correlating aqueous vapor pressures with the composition of various mixtures 


of H5PO0,, when the concentration of HPO, in the solution is expressed as 


weight %, (i.e., mass (HPO, ) x 100/mass (HPO + nH0)). For a specified 


PO, to the 


4 
aqueous vapor pressure, the ratio of the molecular weight of H 


atomic weight of P is 3.16: 


3 


peg loose... 316 
: a oe 
HPO, (wt %) H5PO, wt % 


Taking, as an example, the burning of 1 gram of P in the air at 25°C and 
containing moisture at a partial pressure of 22.39 mmHg, the oxidation to 
phosphoric oxide and hydrolysis of phosphoric oxide will result in 1 x 98/31 
gram or 3.16 gram H5PO,. H5PO, will then hydrate until the aqueous vapor 
pressure of the diluted HAP0, (HPO, + nH,0) equals the ambient vapor 
pressure. At 22.39 mmHg and 25°C, the equilibrium composition of the mixture 
is found to be 20.07 weight % H,PO0,. The YF is then 316/20.07 = 15.7, (i.e., 
under these conditions, the mass of smoke is 15.7 gram for every gram of P 
burned). An additional correction will be necessary if the munition 


efficiency is not 100%. 


DPG (1978) notes that the theoretical curve for phosphorus assumes that 
orthophosphoric acid is produced immediately during combustion. If this 
assumption does not hold, the yield factor may overestimate the effective 
dosage. DPG (1977) also presents a derivation which proves that YF’s are 


insensitive to temperature up to 100". 


Figure B-3 (DPG, 1978) presents yield factor curves based on theoretical 
calculations for P as well as empirical data for red phosphorus and theoretical 
calculations for zinc oxide-hexachloroethane-aluminum (HC). It can be seen 
from these figures that the slope of the lines for all three materials, white 
phosphorus, red phosphorus, and the zinc based munition, is relatively shallow 
over the range of relative humidity from approximately 20% to 65%. For 
example, over the range of 20% to 65% RH, the yield factor for red phosphorus 
increases by a factor of 1.4. However, the theoretical curves rise sharply at 
humidities greater than 65%. Between humidity levels of 65% to 80% red and 
white phosphorus increase by a factor of aoproximately 1.3 and hydrated zinc 
chloride by a factor of 1.7. For field tests conducted on one day, the 
relative humidity would not be expected to vary over this wide a range, so that 


Yield Factor (units) 
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Curve A: Yield factors for the conversion of P to hydrated H;P0, (ref. 8). 
Curve B: Yield factors empirically detenained for RF, 
Curve C: Yield factors for the conversion of Zn to hydrated ZnCl. (ref. 5}. 


Figure B-3. Yield factors as a function of relative humidity for various smoke 
producing agents. Source: DPG, 1978. 


the yield factor between burns should not vary substantially. For example, an 
increase from 40% to S0% RH only produces a difference in the yield factor of 
1.04 for red phosphorus. During the field experiments, which took place on 
several different days, the relative humidity varied from 40% to 80%, but 
within a given day the RH only varied by 2%-4%. 


3.2 Approach to Estimating Uncertainty in the Emission Rate 


A direct method for determining uncertainty values would be to compare 
modeled emission rates with actual emission rates; however, emission rate 
measurements from the field experiments are not available. In the absence of 
actual data, uncertainty values can be modeled based on estimated uncertainty 
values for the input parameters and on the uncertainty inherent in the 
emissions rate model. To estimate uncertainty values for the variables, there 
should be sufficient data available to calculate the probability density 
functions (i.e., measures of variability) of the input variables. To model the 
uncertainty in the emission rate model, an assumption must be made regarding 
the independence or dependence of the variables in the model, as there are 


different approaches for each case. 


The following discussion presents both a qualitative and quantitative 
assessment of the emission rate uncertainty of the various smoke munitions. 
The qualitative assessment highlights the sources of uncertainty of the 
input parameters to the emission rate model and the limitations imposed by the 
amount of data (for each type of munition) available for the analysis. The 
quantitative assessment presents a model for estimating the emission rate 


uncertainty and an example application for several munitions. 


Two factors limit the uncertainty analysis of the smoke munitions data. 
First, to few experiments were conducted for each type of munition to 
determine robust probability density functions for the various parameters. For 
both wind tunnel and field experiments only two burns were conducted for each 
type of submunition, thus means or variances will not be stable or robust. 
Confidence intervals are large when n = 2. Therefore the quantitative 
uncertainty analysis presented subsequently is neither a precise nor an 


accurate measure of uncertainty, but the method is valid with sufficient data. 
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Second, the data from the wind tunnel experiment are suspect. For five of 
the nine munitions tested in the wind tunnel, the mass of aerosolized active 
ingredient, Me as measured by the bubblers exceeded the mass burned, AM, as 
measured by the load cell on which the munition was placed. This clearly 
violates the law of conservation of mass. There are two plausible hypotheses 
for these results. First, the design of the wind tunnel is suspect. [t has 
been observed that there may have been insufficient downwind distance from the 
load cell to establish laminar flow at the sampling point (Carter, (1991), 
personal communication). Second, the wind speed was measured at only one 
location and the assumption of constant velocity may not be correct (Bowers, 
(1991), personal communication). Either or both of these situations could have 
led to an overestimation of air velocity at the point of sampling. Although 
this condition most likely existed for all tests, it was not evident in the 
data from the three types of zinc oxide-hexachloroethane-aluminum (HC) smoke 
munitions (155mm M1 and M2 and the 10Smm canisters) and one red phosphorus (RP) 
munition (81mm Navy wedge). (This observation is based on data reported in 
DPG, 1978.) The reason for the differences between munitions is not readily 
evident, and it is not possible to ascertain given the age of the data. The 
error associated with the inaccuracy of the bubbler data will affect all field 
experiment data because of the fundamental assumption in the emission rate 
model. Namely, the MYF relationship (the ratio of the mass of active 
ingredient aerosolized, ML to the initial mass of the munition, Ms). as 
measured in the wind tunnel experiments, is assumed to apply to the field 
experiments. As the wind tunnel tests overestimated the MYF, the emission 


rate, Q,. would also be overestimated. 


A simple analysis was performed to (1) determine if there was a consistent 
bias in the data which could be used to adjust the mass aerosolized, Me and 
the MYF for those munitions where My. was greater than the mass burned, AM, and 
(2) evaluate the within munition variability for M and AM. The computation 


performed was: ; ‘ 
bias =1- M OM 
where a negative value indicates the mass of Zn or P aerosolized is greater than 


the mass burned and a positive value indicates this relationship is reversed. 
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‘The data were taken from DPG (1978) and the results are shown in 

Table B-1. The first five entries in the table are for the munitions where ML 
the mass of Zn or P aerosolized, exceeds AM mass burned; the next four entries 
are for munitions where M. is less than 4M. For one munition, the 155mm RP 
Navy Wedge, in one test result My was less than AM and for the second test 

MY Was greater than AM. The reason for this anomaly within a munition type is 
not clearly evident. It can be observed from this table that for munitions 
where the mass of Zn or P is less than the mass burned, the within munition 
variability is relatively small, less than 20%. However, for the munitions 
where the mass of Zn or P is greater than the mass burned, the within 
variability was much larger, ranging from approximately 45%-88%. Because the 
within munition variability is so large for those cases when MS is greater than 


AM, it is not appropriate to adjust MY or MYF for these munitions. 


Three conclusions can be made based on this analysis. First, those 
munitions for which the mass of Zn or P aerosolized exceeds the mass burned 
will be dropped from further analysis because the data violate the law of 
conservation of mass and the large within munition variability precludes 
adjusting the original data. Second, the results indicate that there may be 
another source of error affecting the data in addition to the wind tunnel 
design. One would anticipate that the design error would affect all munitions 
equally, however, this is clearly not the case. There could appear to be other 
sources of variation related to the type of munition (e.g., mass aerosolized 
is less than mass burned is for all of the zinc canister munitions). Or there 
may have been some other source of experimental error associated with the 
operation of the wind tunnel-or in the analyses of the bubblers. These 
additional sources of error cannot be evaluated further due to the age of the 
data. Third, if one assumes that the affect of the wind tunnel design on the 
munitions 1s constant across all experiments, then it must be concluded that 
MY. the mass of 2n or P aerosolized, is overestimated for all munitions even if 
MY is less than AM, the mass burned. [It follows then that the MYF, which is 


the ratio of mass aerosolized to mass burned, is overestimated for all 
munitions. 


To estimate the uncertainty in the emission rate model the issue of 


independence or dependence of the variables must be considered. The 


. 
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2.75" WP Wedge 
6" WP Wick 

3" WP Wick 
155mm RP 


Navy Wedge 


81mm RP 
German Wedge 


15Smm HC Mi 
Canister 


155mm HC M2 
Canister 


105mm HC 
Canister 


81mm Navy RP 
Wedge 


Zn or P 
2 M, 


Ag 
A10 


AS 
A6 


A7 
A8 


A15 
Al6 


Bias = 1 - in 


Table B-1 
Evaluating Wind Tunnel Bias 


1 


Initial Mass burned Mass 
Mass (4M) Aerosolized 
(M,) (mM) 

217 125 150 
203 116 129 
102 59 64 
101 §1 77 
63 40 45 
52 36 Void 
114 78 72 
117 69 79 
30 15 24 
35 14 15 

3450 2330 425 

3572 2412 382 
2055 1162 369 
1987 1182 391 
1178 646 154 
1177 633 234 

117 73 64 
116 68 61 
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assumption of independence allows one to use a relatively simple expression 
to calculate overall model uncertainty as the covariance terms can be 
Neglected. However, if this assumption is incorrect and the variables are 


dependent, the estimated uncertainty will be unrealistically low. 


In the case where the variables can be assumed to be independent, 
Goodman (1960) presents the expression: 


N 


—2.2 -2 2 —,2 
o =x vy ae a (xy) [e% + ov] : (B-2) 


where a, = variance of the product xy 


xX  ™ mean of x 


2 Lx - x)? 
c = variance of X = ————. 
x n-i 


y = mean of Y 


a = variance of Y 


2 
Cc 
cv? = squared coefficient of variation of x, oes 


“i 


‘2 
Cc 
cv? = squared coefficient of variation of y, a 


“<I 


The variables in the emission rate model (see equation B-1) are not 
independent, both the munition yield fraction, MYF, and the yield fraction, YF, 
depend on the type of active ingredient. The variance of the product of two or 
more variables, which are not independent, is given as (Goodman, 1960): 


2 


2 
of aX ors y + 2x y(COV,.) 


2 
xy y 
where the variables are defined as above and COV ey? the covariance of xy, is 


expressed as: 


(B-3) 
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To apply Equation (B-2) to the model for Q,, the emission rate, integrating over 
the total burn time, it is written as: 


b 
Q = dtQ, = M,+MYF+YF (B-4) 
° 


Let E = MYF and Y = YF. Assuming the uncertainty associated with the initial 


mass, M is small (i.e., the load cell is accurate), Q, can be divided by M 
and expressed as q 


= EY . (B-5) 


+ 2E Y COV (B-6) 


where COVEY is the covariancé of the munition yield fraction, now expressed as 
E, with the yield factor, now expressed as Y. Since the yield factor is a 
function of the relative humidity (Y = Y(h) where h is fractional relative 
humidity), this value can be incorporated into the expression for 0 as 
follows: 


ee de ae Se 
oq = ¥Y oe + EY o.* 2E YY COVER (B-7 ) 
The fractional relative humidity, h, is = 0 at zero humidity and 1 at 100% 


humidity and Y’ is defined as: 
Y’ = ay (B-8) 


As there are only two values for each munition, Y’ is simply the slope of the 
straight line between the two points. The term cove. is the covariance of the 
munition yield fraction (E) with relative humidity (h).- 
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3.3 Results of Uncertainty Analysis 


Table B-2 presents a sample set of data taken from DPG (1978) and used to 


calculate the 


variance of the quantity QM. expressed as on. The munitions 


included in this analysis were only those where the mass of Zn or P 


aerosolized, ML was less than mass burned, AM, as discussed previously. An 


example of how the variance of QM, (or on) is calculated is given below for 
the zinc munition 155mm HC M1 canister. 


Step 1. 


COVE 


COVE 


Step 2. 


Calculate the covariance of the munition yield fraction (E) and 
relative_humidity (h) for trials B1Ri and B2R1. The values for E 
and h are. given in Table B-2. Equation B-3 is the definitional 


formula for covariance; the computation formula is: 


_ N(Z E+h) - (SE) (Zh) 


n2 


= 2((0.12 x 0.76] + [0.11 x 0.721) - ({0.12 + 0.111[0.76 + 0.721) 
; 32 

= 0.3408 - 0.3404 
— 


1x 1074 


Calculate Y’ according to Equation B-8 where y, the yield factors 
(YF), are 5.2 and 5.0, and h, the fractional relative humidity, is 
0.76 and 0.72. 
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Table B-2? 


Variance of the Ratio of the Emission Rate to the Initial 
Mass, (Q,™M,) for Selected Munitions 


c 
Munition Relative Yield Variance of 
Munition Trial Yield Fraction Humidity Factor Qt 
we 
‘ ___ eeeeeeeeeSeSsSeSesSa ee 
185mm HC M1 . BRI 12 76 5.2 -3 
is Canister B2R1 li 72 5.0 2.15 x 10 
155mm HC M2 B3R1 18 75 5.2 -3 
Canister B4R1 .20 75 5.2 5.41 x 10 
10Smm HC B15 .13 73 5.1 -2 
Canister B16 .19 73 - 5.1 4.68 x 10 
81mm Navy RP Bil .55 70 6.6 
Wedge B12 .52 80 7.8 0. 160 


a. Data taken from Tables III and IV in DPG (1978) for the field experiments 
conducted at the horizontal grid. 


b. MYF values are from the wind tunnel data for the particular munitions. 


c. Estimated values from Figure B~3 for 2nCl, (Curve C) and for the theoretical 
curve for white phosphorus and red phosphorus (Curve A). 
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Step 3. Calculate o according to Equation B-7. The variance of E, the 
munition yield fraction, is 


o =2x 1074, Equation B-7 is as follows: 


and by making the following substitutions: 


o = (8.1)7(5 x 107°) + (0.115)7(5)*(8 x 1074) + 2(0.115) (5.1) 
(5)(1 x 1074) 
=2.15x 10° 


Interpretation of the variance estimates is limited by the potential error 
in the munition yield fraction (MYF) values and the limited amount of data (n = 
2) for each type of munition. However, some summary comments can be made about 
these data and this approach to calculating variance. The variance of the 81mm 
Navy red phosphorus wedge is larger than the three zinc based munitions (zinc 
oxide-hexachloroethane-aluminum). The cause of this difference could be due to 
the type of munition. For example, the red phosphorus MYF is larger than the 
MYF for zinc munitions. Further, for phosphorus munitions the relative humidity 
on the two days of red phosphorus tests differed by 10% where as for two of the 
zinc munition tests it was constant. If there is no variation in relative 
humidity the variance of q is the product of the mean of the yield factor and 
the variance of the MYF (see Equation B-7). 


Bevington (1969) defines the standard deviation, o, as the estimated error 
or uncertainty of a parameter. When a parameter, e.g., x, is a function of two 
or more other variables (x = ab) the standard deviations of a and b, when 
combined, give the uncertainty of x. To put the standard deviation values for 
each munition in perspective, Table B-3 presents the mean emission rates 
integrated over time, the standard deviation, ¢, (as calculated according to 
Equation B-7), and the standard deviation expressed as a percentage of the mean. 
This table shows that the emission rates of the two 155mm HC (zinc) canisters 
(M1 and M2) have the least uncertainty and the 105mm HC (zinc) wedge has the 
greatest uncertainty. The uncertainty of the 8imm Navy RP wedge is only 
slightly greater than the 155mm HC canisters. 
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Table B-3 


Comparison of the Mean and Standard Deviation 


M 


Q 
of the Ratio of the Emission Rate to the Initial Mass | 
° 


Qa Standard Deviation 
Munition Mean — Q Standard Deviation 
M t ‘ 
° of a as a % of the mean 
° 
155mm HC M1 0.587 0.046 7.8 
Canister 
155mm HC M2 0. 962 0.074 7.7 
Canister 
105mm HC 0.817 0.216 26.6 
Canister 
8ilmm Navy RP 3.844 0. 40 10.4 
Wedge 
a Q, . 
al Q, is integrated over the burn time, thus the quantity 7 is actually 
° 
oe 
M b 
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An approach to evaluating emission test results from future experiments 


would be to compute the 95% confidence interval (95% CI = x + 1.960) of the 


Q 
quantity ms for each munition. If a subsequent experiment produced a Ris value 
° ““e 
outside of the confidence interval, than it could be concluded that the 
munition was significantly different from the original sample of munitions. 
Ideally, a larger dataset than that currently available would be necessary to 


assure that the standard deviations were stable. 


An interesting comparison can be made between these results and the 
findings presented in Methodology Investigation Final Report, Validation of 
Transport and Dispersion Model for Smoke (DPG, 1979). In DPG, (197°), modeled 
concentrations along the line of sight (CL) were compared to field measurements 
as well as comparisons between measured and modeled concentration line 
integrated dosage (CLID). DPG (1979) reports that the CLID were consistently 
underpredicted for the zinc based munitions (HC) by a factor of 1.4-1.5 and 
ove-nredicted by a factor of 2 for the red and white phosphorus munitions. A 
comparison of these air dispersion modeling results to standard deviations of 
the source term (presented in Table B-3) shows a similar relationship among 
munitions. The zinc based munitions exhibited less variability than the red 

-phosphorus munition. One hypothesis which could partially explain the 
discrepancies in the modeled versus measured comparisons in DPG (1979) is 
uncertainty in source terms. The red phosphorus has the largest uncertainty 
(as measured by the standard deviation) and it also shows the greatest 
difference between measured and modeled values. The large uncertainty in the 
red phosphorus source term may account for a portion of the uncertainty in the 


modeled red phosphorus concentrations. 
4. Conclusions and Recommendations 


The smoke munitions data from DPG (1978) have been reviewed and a method 
for estimating uncertainty, as expressed as the standard deviation, in the 
emission rates has been presented. Available data for the present analysis is 
limited, but the approach to variance estimation is applicable to larger data 
sets. The quality of the smoke munitions data is suspect due to design 
problems with the wind tunnel. The number of tests conducted on each munition, 


n = 2, means the variance computations are not stable. The uncertainty 
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Oe 


analysis of the source terms indicated that the 105mm HC munition has the largest 
uncertainty, being slightly more than 3 times the uncertainty of the 155mm HC 
munitions and 2.6 times the 81mm Navy RP wedge. For the two 155mm HC canisters, 
the relative humidity was constant for the M2 test runs and varied for the Ml 
test runs. However, when uncertainty is expressed by the standard deviation as 


a percent of the mean, the two munitions have almost the same uncertainty. 


Two recommendations can be made based on the preceding analysis. 
First, the munition yield factor should be determined in an accurate and 
precise manner. This variable is critical, as it is applied to all the field 
test data. Second, to assess meaningfully the uncertainty in the emission 
rate values, more than two sets of data must be available for the results to 


be more stable, and thus reliable. 


To use this technique of computing variance to evaluate munition 
variability, the following plan is suggested. Assuming the wind tunnel design 
issue has been resolved and accurate MYF’s determined, field data from numerous 
test burns should be collected. This should be done at the horizontal grid, 
but without simultaneous aerosol and meteorological monitoring. This is 
because to compute uncertainty for the emission rate value, the only parameters 
needed are initial weight, Mo. munition yield fraction, and yield factor (as 
computed based on relative humidity); thus only Mo and “%‘RH would need to be 
recorded. Numerous burns for each type of munition should be conducted, e.g., 
N = 30. Tests should be conducted under as wide a range of humidities as 
possible. Such a database could then be used to compute robust variance and 
95% confidence intervals. Subsequent test results would be compared to this 
confidence interval to determine if the emission rate is statistically 
different. Further, with sufficient data to compute stable variances, these 
values could be incorporated into an expression with meteorological variances 


to compute overall model uncertainty. 
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Appendix B-2 
Analysis of Fog~Oil Smoke Emissions 


1. Introduction 


This analysis focuses on the fog-oil smoke generator emission rates as 
tested at the Dugway Proving Ground in March and April, 1985 (Liljegren et 
al., 1988). Liljegren et al. (1988) reported the details of how the 
experiment was conducted. Only the method of computing the emission rate will 
be addressed here. The objectives of this analysis are threefold: (1) to 
assess the variability in the emission rate from the fog-oil smoke generator, 
(2) to compare two methods of computing the emission rate, and (3) to evaluate 


the influence of averaging time on emission rate variability. 
2. Method 


The fog-oil smoke emission data were taken from Liljegren et al. (1988). 
The configuration of the smoke test provided for the "instantaneous" 
measurement of the emission rate based on the weight loss of the oil drum and 
the exit velocity of the generator. Due to mechanical difficulties, the 
weight loss values are actually i-minute averages, as opposed to 
“instantaneous” values. Plots of these data (and exit temperature) versus 
time for seven experiments were reported in Figures 3.2 through 3.8 in 
Liljegren et al. (1988). (For reference, copies of those figures are 
included, with the original figure numbers.) These data were digitized and 
the means, x, and standard deviations, o, were computed for each experiment. 
These values and the coefficient of variation (¢/x) for each test are 
presented in Table B-1. Some experiments were screened to remove values which 
were approximately zero at the end of the test. Otherwise, the data reflect 


the actual emissions for the duration of the experiment. 


Liljegren et al. (1988) computed a time integrated emission rate based on 
the total mass of oil burned divided by the duration of operation of the smoke 


generator. These values, taken from Table 3.1, are shown in Table B-2. 


To evaluate the influence of averaging time on the variability in the 
emission rate, the emission rates were re-calculated based on averaging times 


up to 510 seconds. To do this the data were processed through a program which 
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Average Instantaneous Standard Coefficient 
Figure No. Test No. Emission Rate (g/s) Deviation of Variation 
3.2 TOOO3 21.7 10.9 0.5 
3.3 T0004 33.8 8.9 0.3 
3.4 Tooos 21.5 5.9 0.3 
3.5 TOO06 31.1 4.2 0.1 
3.6 T0007 30.9 7.5 0.2 
3.7 TOOOs 37.9 7.5 0.2 
3.83 TOOO9 37.9 16.8 0.4 
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Table B-1 


Variation in Instantaneous (i.e., one-minute averaged) 
Fog-Oil Smoke Emission Rates 


Table B-2 


Comparison of Integrated and Averaged Instantaneous Methods 
of Calculating Emission Rates 


Integrated « Average Instantaneous Percent 
Fisure No. Test No. Emission Rate (g/s) Emission Rate (g/s) Difference 


3.2 T0003 24.3 21.7 -12 
3.3 70004 22.5 33.8 +33 
3.4 T0105 22.4 21.5 -4 
3.5 To0V6 28.7 31.1 +8 
3.6 To00C7 30.7 30.9 -0.1 
a, TOOC8 36.3 37.9 +6 
3.8 Toco9 43.2 37.9 =12 


* . 
Taken freu Liljegren et al., (1988), Table 3.1 
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interpolated the emission rate at 30 second intervals from the existing data 


which are at approximately one minute intervals. From this processed data, 
means and standard deviations were computed by grouping the data at 30 sec, 60 
sec, 90 sec, etc. increments up to 510 seconds (8.5 minutes). To evaluate how 


the variability changed with averaging time, the following equation was used 


7 sec 
R’ = [ oe (B-1) 
30 sec 
where Tien is the standard deviation for the i-th averaging time (i = 60 sec, 


90 sec, 120 sec, 510 sec) and co is the standard deviation from the 


30sec 
30 second averaging time. 


3. Results 


The first objective of this analysis is to evaluate the variance in the 
fog-oil generator emissions. Table B-1 shows that the variation in the 
emission rate, expressed as the standard deviation, can be large. The 
coefficient of variation (CV) normalizes the standard deviation by ‘the mean. 
Test TOOO06 has the smallest CV - 0.1. Tests TOOO3 and TO009 have the largest 
CV’s, 0.5 and 0.4 respectively. During these tests the smoke generator 
malfunctioned and the emissions were particularly erratic (see Figures 3.2 
and 3.8). 


The secord objective of this analysis is to compare the methods for 
calculating the average emission rate. Table B-2 shows a comparison of the two 
methods of computing the mean emission rate. Liljegren et al. (1988) 
calculated the time integrated average by dividing the mass of oil burned by 
the duration of operation of the smoke generator. Our approach averaged the 
“instantaneous” 1 minute emission rates, as digitized from plots. The percent 
difference between the two methods ranges from -12% to +33%. The largest 
difference in the computed emission rates is for experiment T0004. A review 
of Figure 3.3 and the digitizing indicatés that there were no errors in that 
step and that the mean rate of 33.8 g/sec accurately reflects the data in the 
plot. [It is hypothesized that there may be an error in the field data 
relating to oil burned or duration of operation of the smoke generator, which 
may have caused Liljergren et al. to compute a lower emission rate. (In 
addition, Liljergren et al. note that the data for this experiment were 


divided into two groups based on a shift in wind direction. [t is not clear 
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how this was done of if this would affect the emission rate calculation. The 
smallest difference was for TOOQ07, an experiment when the smoke generator 
functioned in a very consistent mode (see Figure 3.6). There does not appear 
to be a bias in the Liljegren et al. method, with their approach producing a 
higher mean rate four times and a lower rate three times. The r? for the two 
methods of computing the emission rate is 0.56. This means the time 
integrated method of computing emissions can only account for 56% of the 


variance in the method based on averaging the “instantaneous” values. 


The third objective of this analysis is to evaluate how the variablity in 
the emission rate changes with averaging time. The data from five 
experiments, TOOOS (Figure 3.4) through TOOO9 (Figure 3.8), were analyzed 
using equation (B-1). Two experiments, TOOQO3 and T0004 (Figures 3.2 and 3.3, 
respectively), were excluded from this analysis due to insufficient data. 
Table B-3 presents the squared standard deviation ratios by averaging time for 
the five experiments as well as the median ratio for each averaging interval. 
These data are plotted in Figure B-1, with the solid line representing the 


median value and the minimum and maximums represented by the coxes. 


A theoretical curve is also plotted in Figure B-1, representing the 


equation: 


2 
o (T.) 1+ 30,/2T, 


- TST (B-2) 
a I 


o* (30s) 


where T; is the integral time scale of the physical process. This equation is 
based on work by G.I. Taylor in the 1920’s and is further discussed in Section 
VIII of the main text. This theory assumes that the time scales of the 
fluctuations cover a wide range of values, and the time series data in the 
figures verify that this is the case. An integral scale, qT of 150 seconds 


provides the best fit to the data in Table B-1. 
4. Conclusions 

This analysis shows that the variation in the one-minute averaged 
emission rate from the fog-oil smoke generator can vary from 10% to 50% of the 
mean rate, based on the coefficient of variation. From the comparison of the 


two methods of computing the emission rate, it appears that when the smoke 
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Variations in Emission Rate with Averaging Time 
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Figure B-1. Variation of emission rate fluctuation variance with averaging 


time. The solid line represents the median of the five 
experiments and the stars represent the range at that averaging 


time. The dashed line represents the theoretical curve given by 


Equation (B=-2). 
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generator is operating in a consistent manner, both computational methods 
produce the same result (e.g., TOOO7, Figure 3.6). However, when the smoke 
generator malfunctions (e.g., TOQ0Q3 and T0009, Figures 3.2 and 3.8, 
respectively), emissions can be highly variable, and the method of Liljegren 
et al. will over estimate emissions relative to our method of averaging 
“instantaneous” emissions. An analysis of the emission rate variability with 
different averaging time showed that a time interval of 150 seconds best 


represented the experimental data. 
5. References 
Liljegren, J.C., W.E. Dunn, G.E. DeVaull, A.J. Policastro, Field Study of 


Fog~-Oil Smokes, Supported by U.S. Army Medical Research and Development 
Command, Fort Detrick, MD, January, 1988. 


Intentionally Blank 


APPENDIX C 


UNCERTAINTIES IN SOURCE EMISSION RATE ESTIMATES 
USING DISPERSION MODELS 


Reprint from Atmospheric Environment, 24A 2971-2980 (1990) 


 Alewsphene Exeveument Yok ba. No tl pp 97) rah ie 
Printed 1a Great Baas, 


UNCERTAINTIES IN SOURCE EMISSION RATE ESTIMATES 
USING DISPERSION MODELS 
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Abstract—The source emission rates during the Prairie Grass dispersion experiments were careiully 
observed and were adjusted by the expenmentalists so that they were about twice as high during unstable 
conditions as dumng stable conditions. The question was asked whether observed concentrations and 
meteorological conditions could be used in dispersion models in order to predict source emission rates and 
veniy this factor of two difference. Three cypes of simple dispersion models were applied to this problem. 
with the resuit that for the model based on Monin—Obukhov simuilanty theory, the uncertainues in 
predictions of source emission rates for individual runs were at best about = 10-20% when voserved cross- 
wind integrated concentrations from the 50m arc were used. Consequently this model could discern the 
factor of two difference in average source emission rates for the two sets of fieid tmals which consisted of 
about 20 runs each. However, some models. such as the Gaussian plume model. exhibit uncertainties of 
about = 70% to a factor of two in predictions for individual runs. and hence could not discern the difference 
in average source emission rates when concentration observations at downwind distances of 100-300 m are 
used, It is found that the use of observed cross-wind integrated concentrations produces more accurate 
conclusions than che use of observed point concencrations, for the uncertainties in predictions of source 
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emission rates are about a factor of two larger when the observed point concentrations are used. 


Key word index: Dispersion models, uncertainues in models. source emission estimates. 


OBJECTIVE AND METHODS 


AS source emission rates for air pollutants are seldom 
weil-known, air pollution control decisions must often 
be made based on observations of air pollution con- 
centrations on monitoring networks. Observed con- 
centrations can be combined with meteorological 
observations (e.g. wind speed and direction, stability, 
and mixing depth) and used in so-called hybrid 
source-receptor models in order to estimate source 
emission rates (Watson, 1989). These models are based 
on the concept that atmospheric transport and disper- 
sion modeis are a mathematical link between observed 
concentrations and predicted source emission rates. 
The specific objective of this research, funded by the 
U.S. Army, is to develop methods for estimating 
whether there is a significant difference in the source 
emission rate of similar types of sources from one 
experiment to another. For example, a typical source 
might be a fog oil generator (used for smoke obscura- 
tion purposes), and the question may be whether the 
source characteristics of the fog oil generator have 
significantly chaaged between 1980 and 1990. It is 
assumed that the source emission rate is not directly 
measured, but that several field trials are carried out at 
the same site in 1980 and in 1990. In each field trial, 
cross-wind integrated concentrations are observed at 
a distance 100 m downwind of the source, and wind 
velocities, turbulence, and vertical temperature gra- 
dients are observed on a tower near the source. The 
source position and orientation is not varied. Some 


. 


data of this type exist, and the procedures iil even- 
tually be tested with these data: however, cirst we have 
chosen to apply the source emission rate estimation 
methods to a simpler data set. 

The Prairie Grass dispersion experiments produced 
a high-quality database (Barad, 1958) that has been 
used extensively in the development and testing of 
dispersion models. The experiments were simple, with 
near-ground-level continuous point sources, and ex- 
tensive concentration and meteorological observa- 
tions were made over the dat. homogeneous terrain. 
Source emission rates were also reported, and were 
deliberately varied so that they were twice as hugo 
during the day as during the night. Because of the care 
with which the Prairie Grass expenments were con- 
ducted, this database serves as an excellent test bed for 
study of the uncertainty in source emission rate es- 
timation procedures. If the procedures cannot discem 
the known factor of two difference in average source 
emission rates in this highly controiled set of ex- 
periments, it would be of little use at more complicated 
sites, where the terrain may be more complex and the 
observations are likely to contain more errors. 

Three representative dispersion models are used to 
relate observed concentranons and meteorological 
conditions to source emission rates. The three models 
include a statistical regression model, a Monin- 
Obukhov similanty model, and a Gaussian plume 
model, The Student-c test is used to determine whether 
the models can discern a factor of two difference, at the 
95% confidence level, between average night and day 
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source emission cates. Another output of this proced- 
ure ig an estimate of the typical error or uncertainty 
the source emussion rate estimate. 


DESCRIPTION OF PRAIRIE GRASS DATABASE 


Barad (1958) discusses the Prairie Grass database in 
great detail. SO, tracer gas was released over periods 
of about 10 min from a point source located at an 
elevation of 0.45 m. Ten minute averaged tracer con- 
centrations were obtained from measurements ob- 
served at an elevation of 1.5m along arcs at distances 
of 50, 100, 200. 400 and 800m from the source. 
Supporting meteorological observations were made 
from a nearby tower. tocated in a fat area, represent- 
ative of the site, where the average surface roughness. 
2, was determined to be 0.6cm. Data from 44 ex- 


periments were used in our analysis. covering a wide. 


range of meteorological conditions. The data are 
approximately equally divided into unstable and 
stable periods. 

The Prairie Grass data are listed in Table 1, based 
on information in the report by Barad (1958) and in 
papers by van Ulden (1978), Nieuwszadt (1980), and 
Briggs (1982). The 2 m wind speed. u. the 2 m standard 
deviation of wind direction Auctuauons o,, and the 
16m to 2m temperature difference. OT. were given in 
Barad’s i198) report, and the inction veiocity, u,, and 
the Monin—Obukhov length, L, were esumated by 
van Ulden (1978). The mixing depth and convective 
scaling velocity, w,, were calculated for most unstable 
experiments by Nieuwstadt (1980). The stability ciass, 
SC, is estimated using observations of surface rough- 
ness, z,, Monin—Obukhov length. £, and Goider’s 
(1972) nomogram. Briggs (1982) carefully reviewed all 
of the dara and suggested that certain arcs for some 
experiments be removed from the data set because of 
problems with those particular data. At each down- 
wind distancs, the table lists observations of C/Q, 
C”’Q, and o,, where C is the maximum concentration 
on that arc, C” is the cross-wind integrated concentr3- 
tion and a, is the standard deviation. of the lateral 
concentration distribution. The C’ and a, data were 
given by Nieuwstadt (1980) for most of the unstable 
periods. The standard deviation of vertical wind direc- 
tion fluctuations, ¢,, has been estimated using bound- 
ary layer similarity theory. 


DISPERSION MODELS 


Data from the Praine Grass experiments i Barad, 
1958) were used extensively by Pasquill (1961) and 
others in the development and testing of a Gaussian 
diffusion mode! now known as the Pasquill—Gifford— 
Turner model (Gifford, 1962, 1968, 1976: Turner, 
1967), which is the basis for most decisions regarding 
air pollution controi in the U.S. The data were also 
part of the database used by Nou (1963) in the 
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development of the empirical OB, OG model, which is 
the basis of the U.S. Air Force procedures for calcu- 
lating toxic gas impacts. Because the zariy data ana- 
lvses did not make use of Mflonin—Obukhov simiuanty 
modeling or convective scaling conceots, there was 2 
flurry of activity with reanaivses of the Praime Grass 
data from this new point of view in the late 1970s 
(e.g. van Ulden, 1978; Horst. 1979: Nieuwstadt, 1980: 
Venkatram, 1981: Briggs, 1982). In several of these 
Papers, the new scaling parameters (e.g. the mixing 
depth, A, the friction velocity, u,, and the 
Monin-Obukhov length, L) were desived by reanaly- 
zing the original field dara. 

The variety of dispersion models in the references 
listed above can be grouped into three classes. 


Class |. Empinecai or statisticai regression moceis 
(e.g. Nou, 1963). 

Class 2. Similarity models (e.g. Briggs, 1982). 

Class 3. Gaussian piume modeis (e.g. Gifford. 
1976). 


Because Nou’s (1963) empirical or statistical regres- 
sion formula was based on several other databases 
besides the Prairie Grass database, and he did aot 
suggest a formula for the cross-wind integrated con- 
centration, we decided not to use his formula directly, 
dut instead we applied a muit:vamate linear regression 
procedure to the data in Tabie i :n order to demve tne 
foilowing best-it sower-iaw formuias: 


C/Q = 0.000137x ~'-$1(DT— 10°F)*:*? () 
C/Q =» 0.00666x ~ '-°3(DT+ 10°F)?-$* (2) 


where C/Q and C’/Q are in units of sm~* ands m7~?, 
respectively, x is in units of m and D7 is in units of °F. 
This is the same statistical procedure applied by Nou 
(1963), and the units for all variables are consistent 
with those that he used. Nou (1963) aiso used the 10°F 
additive factor, which is necessary to prevent negauve 
values of the temperature term. These formulas ex- 
plain 91 and 84°%, respectively, of the vamance in the 
C,Q and C”’Q observations, where most of the var- 
ance 1s explained by the x term. Note that these 
formulas are aot based ar all on physical insights, dur 
are based on least-square fits of linear reianionships to 
the data. 

The group of similarity models proposed by van 
Ulden (1978), Horst (1979), Nieuwstadt (1980), Venka- 
tram (1981) and Briggs (1982) is based on appiicauons 
of Monin—Obukhov similarity theory and convective 
scaling similanty theory. According to Monin~ 
Obukhov similanty theory, the following functional 
relations can, be postulated for dispersion from contun- 
uous ground-level point sources in the surface bound- 
ary layer (Bnggs, 1982): 


Cu,x*/Q =f (x/L) (3) 
Cu x/Q =fz(x/L) (4) 


where f, and f, are non-dimensional! universal func- 
uons. The friction velocity, u,, and the Monin- 
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Obukhov length. © are the fundamentat scaling veio- 
city and scaling length in the surtace boundary layer, 
according to Monin—Obukhov similanty theory. The 
Praine Grass data are plotted in Fig. | in the dimen- 
sioniess forms suggested by Equations (2) and (3), 
ilustrating that the data are ordered quite weil by 
these similarity relations. 

Because the data in the four parts of Fig. | appear to 
approach a constant as x/L approaches zero. the 
following functional relation is proposed: 


S(x/L) ma(l + bx/L¥. (5) 


A least-squares algorithm was applied to each figure, 
with the following results: 


A((x/L) = 3.371 —O.019x/L)7 "36 <0 (6) 
flr L) 3.011 = 220x/L)°-5” xL>0 (7 
flr L) 2 1.061 -0.021x/2)7'-7* x L <0 (8) 
Sylx/L) = 1.07(1 + 0.10x/L)9-S4 x/L>0. (9) 


These curves are drawn on the figures and appear to 
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provide a good fit ac ail vaiues of x/L. However. they 
have been denved with no requirements that certain 
physicaily-based asymptouc functuonal relanonships 
be satisfied. For example. Briggs (1982) points our that 
f,(x/L) should be proportional to (—x/L)!? in the 
limit of free convecuon (—x/L—~co). Consequently, 
Briggs’ (1982) suggests a third term in Equauon (5), 
which he claims to better account for “convective 
sweep-out” of the plume at large {x/Z|. The similanty 
moded is represented by Equations (6)-{9) in our 
analysis rather than the equations proposed in any of 
the above references because the other equations are 
based on a slightly different subset of the Prairie Grass 
database. 

The Gaussian plume mode!s comprise the ‘third 
class of models in our analysis. The Gaussian model is 
based on the formulas (Gifford. 1976) 
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(lL!) 


C/Q =(tuo,o,)~' exp( =(z, -2.)° 207) 
C7/Q =(2/)!/2 (ua) ' expt —(z, — 2)? 207). 


The source height, :,, and receptor height, :,, equal 


Stable 


Cu,x?/Q 


188.4 


18.3 


x/L 


Fig. |. Praime Grass data plotted in the dimensioniess form suggested by Equations (3) and (4). Best-tit lines of the 
form atl -6x/L)* are drawn on each figure {see Equations (6)-9)). 


Source emissiva rate Slates 


0.45 m and [.5 m. respectively. in the Prairte Grass 
experiments. The wind speed. u. ts assumed (o be that 
measured at a he:ght of 2 m (see Table |). The values 
for lateral and vertical dispersion coerticieats. o, and 
a,, that are subdstututed inco Equations (10) and (11) 
are based on suggestions bv Brggs (1973), who used 
the Prairie Grass experiments, along with many other 
fieid experiments. in their derivation. These ¢, and ¢, 
formulas are given in Table 2. The required stability 
class is found in Table | for each Prairie Grass test 

Note that the empirical and similartty classes of 
dispersion models have been best-fit to the exact same 
data (Table |) that will be used for further testing. The 
Gaussian plume modei has noc been subjected to these 
same best-fit procedures. but has been based on data 
from many expenments. Because of these differences 
in databases used to derive the models. it is expected 
that the Gaussian plume mode! will have the least 
success in any compamsons with held data from the 
Praine Grass database. 

The dispersion models in Equations (3) and (4) and 
(6)+{10) are inverted in order to predict source emis- 
sion rates: 
Q,2C,(C,/Q), (12) 


Q,202,(C7/Q),, (13) 


where subscmpts o and 5 represent observed and 
predicted vamaoles, respectuvely. Because of the 
characterstics of the muitipie iinear regression or 
least-squares statistical procedures. an equation that 
produces a best-fit for CQ or C”/Q may not necessar- 
ily product a best-fit for Q. These statisucal procedures 
attempt to minimize the mean-square error and to 
force the mean of the observed variable to equal the 
mean of the predicted variable. We should therefore 
not be surprised if an equation which produces zero 
mean bias in C,Q is discovered to produce a mean bias 
of 20-30% or more in Q. 


STATISTICAL PROCEDURES 


A primary objective of this research is the develop- 
ment of methods for estimating whether there is 
a significant difference in the source emission rate 
of similar types of sources from one experiment to 
another, as determined dy observations of meteoro- 
logical conditions and of point concentrations or 
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cToss-wind integrated concentrauons. Dispersion mo- 
dels, such as the three presented above, can be used to 
remove the effects on the concentration ooservations 
of varianions in meteoroioical parameters and down- 
wind position of the concentration monitor. 

A listing of the observed source emissions. Q, during 
the Prairie Grass experiment is given in Table |, and 
the data are plotted in Fig. 2 as a function of the 
stability parameter 1/Z (the point for run 46 has been 
excluded, since it was an evening transition period 
when Q was still high although 1/L had just become 
positive). During the Prairie Grass experiment the 
night-time emissions, Q, of tracer gas were controlled 
so that Q averaged 45.2 gs~' with a range from aoout 
38 to 58 gs~! and a standard deviation of 5.91gs7', 
and the dayume emissions were controiled so thar 0 
averaged 98.25 gs~' with a range from about 90 to 
105 gs~' and a standard deviation of +48 gs~'. The 
expenmentaiists deiiderately maintained this fact’ of 
2 difference in day-night Qs so thatthe magmituc sf 
concentrations at the mvaitors would aot vary 
much from day to night. With 19 stable (night) runs 
and 20 unstable (day) runs, the Student-¢ parameter 
(Panofsky and Brier, 1968, p. 63) can be calculated in 
order to determine if the daytime @ is significantly 
different from the night-time 0: 


where subscripts d and n indicate day and agnt. If iti is 
less than 204, then the difference (G,-Q,) is not 
significantly different from zero, at the 95% confi- 
dence level. In fact, using the observed @ and ag 


a3 a1 
1/L (1/m) 
Fig. Observed tracer gas source emussion rates. Q, as a 


funcuon of inverse Monin—Obukhov length, 1/2. for the 
Prarie Grass data (run 46 exciuded). 


Table 2 Formulas recommended by Snggs (1973) for ¢, and ¢., for rural 
conditions 


Stability class a,(m) 


DANQONDS 


0.22x41 =0.0001x)7 '/? 
0.16xt1 ~0.0001x)7 1? 
O.1 Ley 1 = 0.0001 x) 7 1? 
0.08.x( 1 + 0.0001.) 7 '/2 
0.0611 ~0,0001x) 7 '/? 
0.04 <1 t ~0.0001x)7 1? 


gm) 


0.20x 

0.12 

0.08.xt! ~0.0002x)~ 7 
0.06x¢1 ~0.0015x37 '/? 
0.03.x11 +0.0003x)7* 
0.016x11 =0.0003x)"' 
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figures quoted above. the calculated Student-c para- 
meter is 30.86 for the difference between the mean 
night-time and daytime emission rates, which implies 
that the difference ts significant at far greater than the 
95% contidence level. 

Because there is such a clear difference in the 
observed day and night tracer gas source emission 
rates. Q, at the Prairie Grass site. it is inceresting to ask 
whether the dispersion models would be able to 
predict this significant difference in the mean emission 
rates, based on observed concentrations and meteoro- 
iogical variables. The three types of models described 
earlier were used to prepare predictions of source 
emission rate, Q, (one set of predictions using C, and 
another set using C?). This calculation of Q, was made 
individuaily using data from each of five monitonng 
arcs (50. 100, 200, 400 and 800 m) for C,, and each of 
three monitoring arcs (50, 200 and 800 m) for C2. For 
each set of predicnons, values of aight-time and 
daytime @, and gg, were calculated. Then Equation 
(14) was used to calculate the Student-¢ parameter. If 
the resulting ¢ value exceeds 2.04, it is concluded that 
the model successfully simulated the difference in 
source emission rates. If ¢ is less than 2.04, then it is 
concluded that the mode! has failed to simulate the 
difference. 

Besides invesugating the day-night difference : 
mean source emussion rate predictions, :t is aiso of 
interest to invesugate the abtiity of the models to 
predict emission rates for ail Prairie Grass runs taken 
as a group. For this purpose, che following per- 
formance measures are calculated (where an overbar 
represents an average over ail the database). 


5 
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Relative bias: (0, -0,),. 


raiQ, we Q,) (Q, im 01/72, %,- 
Fractions of predictions, Q,, within a factor of two of 
observations. Q,. 


Correlation: 


Normalized mean square error: (Q,—Q,)*'0,0,- 


These performance measures are tabulated for each 
mode! and for C, and C2 on each monitoring arc. [n 
addition. the individual cesiduals (Q, —Q,) are plotted 
as a function of downwind distance for each model. 


RESULTS AND CONCLLSIONS 


Predictions of average duytime and night-time source 
emission rates 

The procedures for estimating the daytime and 
night-time average source emussion rates reviewed 
above were applied to the Prairie Grass database, with 
the results given in Tables 3 and 4, for C, and CY, 
respectively. By comparing the numbers in the two 
tables, it is evident that use of the observed cross-wind 
integrated concentration, C2, produced better results 
than use of the observed point concentration. C,. The 
standard deviations, oy, of the predicted source em:s- 
sion rates are adout 50° iarger dunng the might anc 3 
factor of 2 or 3 larger dumng the day, for C, than ‘or 
C2. Consequently the calculated Student-: parametes 
are about twice as large for C? than for C,, implying 
that there is more confidence in the conc: isions for C?, 
We expected that there would be a difference, because 


Table 3. Predictions of source emission rate. Q,(gs~'), for night-cume and daytime Praine Grass 

runs, using odserved maximum concentrauons, C,, on momsonng arcs at distances of 50. !00. 

200, 400 and 800 m. Three different models are used to caicuiate @, =C,(C.Q),. The average. @,, 

and standard deviation, ¢9,, for night-time and daytume conditions is listed. Student-r :s 
calculated using Equauon (14) 


Monitonng Night (N = 19) Day (N #20) 
distance Student Q, a9 om Tay 
Model (m) t (gs~') (gs~') igs~') igs7‘) 
Observed Q 50.86 45.2 5.9 98.2 3.5 
Regression 50 7.40 30.1 12.0 144.6 64.7 
100 $.47 $3.8 19.4 139.1 71.6 
200 4,27 52.0 23.0 128.3 725 
400 220 612 27.3 93.4 $6.2 
800 —h19 84.2 $1.6 64.5 48.9 
Similarity sO 5.68 39.9 141 70.8 18.5 
100 4.62 49.2 16.2 34.0 1 
200 5.51 48.7 16.2 A0t1 37.1 
400 §.52 45.3 15.6 105.4 43.2 
300 4.30 46.3 14.6 116.2 60.1 
Gaussian 50 2.40 $5.2 45.8 "15.5 96.9 
100 1.07 76.3 94.9 109.2 92.2 
200 0.05 98.1 180.9 100.2 70.1 
300 —0.88 - 118.8 196.3 76.8 62.7 
300 - 1.67 150.1 264.1 48.0 30.1 
C-6 


Source emission rate Somates 


Table $. Predictions of source emission raic. Q,igs ~') for night-time and dayume Praime Grass 

runs. using OOserved cross-wind integrated concentrauons. Cl, on monitonng ares at distances of 

50. 200 and 300 m. Three different models are used to calculate Q, = C1.(C7'Q),. The average. d,. 

and standard deviation, 9... Jor nignt-cme and daytme conditions listed. Student-< is caicuiated 
using Equauon (14) 


Monitonng Night (N = 19) Day tN =20) 

distance Student @, Say @. Say 

Medei (m) t (gs~') igs~') (gs7') (gs~') 
Observed Q aay by A 45.2 5.9 97.8 4.7 
Regression 50 179 29.7 9.3 140.1 24.8 
200 7.43 46.9 16.9 124.8 40.1 
800 —0.58 69.4 34.3 621 36.4 
Similarity Bit] 16.24 39.7 3.3 89.6 9.0 
“a 9.54 49.0 9.0 100.1 20.3 
300 4.58 45.2 13.3 39.6 50.2 
Gaussian 30 6.18 61.7 20.9 161.0 38.4 
200 2.33 65.5 60.3 116.6 33.9 
300 -0.99 $5.5 145.1 6.9 20.3 


a, is another complicating factor included in C, thac 
would act to increase the uncertainty, but did not 
expect that the difference would be this large. 

The ability of the models to arrive at the proper 
answer (i.e. that the daytime @ is significantly larger 
ihan the night-time @, with at least 95°% confidence) 
can be seen by idenufying cases where © 22.04 in 
Tables 3 and 4. [t is seen that most models are able to 
reproduce this conclusion at most arc distances. How- 
ever, faise conclusions (¢<294) are reachea for the 
following irc and models. 


C,: 100 mare Gaussian 

200 marc: Gaussian 

400 m arc: Gaussian 

800 m arc: Regression and Gaussian. 
C2: 800 mare: Regression and Gaussian. 


It is seen that the use of C, and C? from the closest 
arc (50 m) leads to the correct conclusion for ail three 
models, with a value for Student-c ranging from 2.4 to 
7.4 for C, and 5.2 to 17.2 for C2. The regression and 
similarity models are seen to be better able to simulate 
the difference in observed Q (i.e. they yield a higher ¢) 
than the Gaussian modei. However, the aumber of 
fals= conclusions increases as downwind distance in- 
creases. At the 800m arc in both tables, only the 
similarity model yields the proper conclusion. Figs 3 
and 4 contain extreme examples of plots of predicted 
source emission rates, Q,, as a function of 1/L, in the 
same format as the observations in Fig, 2. The sim- 
larity model predictions in Fig. 3, using observations 
of C2 on the SG m arc. produce patter ; simular to the 
observed data. but with slightly more scatter. How- 
ever, the Gaussian piume model predict:ons in Fig. 4, 
using observations of ©, on the 800 m arc, obviously 
miss the overall trend of the observations of Q, as weil 
as containing much more scatter. 

The failure of the regression modei and the Gaus- 
sian model to estimate @ using C, and C? observations 


Qo (g/s) 


od.k 3.3 aes 3.2 3.3 


L/L (1/m) 


Fig. 3. Tracer gas source emission rates, Q,, predicted by 

the similarity modei as a function of inverse Monin— 

Obukhov length, 1/L, for the Prairie Grass data (run 46 

excluded). Observed cross-wind integrated conceatra- 
tions on the 50m arc are used. 


Q% (g/s) 


3.1 a3 3.1 4.2 3.2 


L/L (1/m) 


Fig. 4. Gas source emissions rate predicted by the Gaus- 

sian plume model, Q,, as a funcnon of inverse Monin- 

Obukhov length. 1/L. for the Praine Grass data irun 46 

exciuc 7d). Observed concentrations on the 300 marc are 
used. 


atthe! est distances may been seen by invesugaung 
the vana..on of the ratio C, C, with distance for the 
three models for stable and unstable conditions (Figs 5 
and 6). The C, C, values for the similarity modei have 
very little trend with distance and are always centered 
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C/Q, Regression 


x (m) 


C/Q, Gaussian 


139.28 
18.28 
Ly al > = 
Tat = 
: +’ P 
3.31 = 
19 31 13aS——s«éiRDSSCi«CLHZ 


x (a) 


Fig. 5. Ratios of predicted to observed concentrations 
(C/Q) for the Prairie Grass dataset, as a function of 
downwind distance, for three models (regression, sim- 
varity, and Gaussian) for stable conditions. V is the 
number of data points at each distance. The midline of 
each box piot is the median, and the other lines 
represent — one and two standard deviauons. 


on 1.0. In contrast, at large distances, the regression 
and Gaussian models tend to underpredict (C, C, < 1) 
during stable conditions (Fig. 3) and overpredict 
(C/C, > l) during unstable conditions. It 1s important 
to note that a model which overpredicts the concen- 
tration, C, will underpredict the source emission rate, 
Q. Consequenuy, at large distances, the regression and 
Gaussian models can be expected to overpredict the 
sOurce emission rate, Q, during stable conditions, and 
underpredict it dumng unstable conditions. But since 
the observed stable @ is about 45 gs~' and the 
unstabie @ ts about 98 gs~*, she Os predicted by these 
two models will tend to be ether nearly the same, or 
their relative magnitudes may even Se switched. This 
problem is seen to occur for the Gaussian model on 
the 800 m arc, using the C, data, where the following 
discrepancy 1s found in Table 3 and can be seen in Figs 
2 and 4 
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Fig. 6. Ratios of predicted to observed concentrauons 
(C/Q) for the Prairie Grass dataset. as a function of 
downwind distance, for three models (regression, sim- 
ilarity, and Gaussian) for unstable condiuons. N is the 
number of data points at ¢acd distance. The mudline of 
each box plot is the median, and che other ‘ines 
represent = one and two standard ceviations. 


Thus the Gaussian mode! leads to the opposite 
conclusion regarding daytime and night-time © differ- 
ences if the 800 m C, data are used. Note that if the 
experimentalists had controlled the emissions in the 
Opposite way, such that the night-time Q exceeded the 
daytime Q, then these trends in the errors in the 
regression and Gaussian models would have led to a 
correct conclusion regarding the differences in Q. 


Night Day 
Observed Q, F 35.2 98.2 
Gaussian predicted Q, 150.1 43.0 


Consider the best-performing mode! in Tabie 6. 
where C} data from the 50m are are used in the 
similanty mode! to predict the average source emus- 
sion rate, Q (see Fig. 3). The caiculated Student-< value 
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Table 5, Evaluations of predictions of source emission rate, Q,(g5~'). for ail Praime Grass runs. using 
observed maximum pot concentrauons, C,, on Monitoring arcs ul distances of 20, 200 and $00 m. Three 
different models are used to catculate Q,=C (CQ), 


Fracuon 


Distance to C, Q-2, Correlation within Q.- Qs 
” observation Mode! d,igs7") 2, r factor of 2 0.0, 
Observed 724 
50m Regression 83.3 0.23 0.77 0.82 0.53 
Similarity $5.7 -0.23 0.74 0.90 0.15 
Gaussian 86.1 0.19 0.36 0.80 0.97 
200m Regression 92 0.26 0.59 0.85 0.51 
Sisularicy 15.6 0.04 0.71 0.92 0.14 
Gaussian 99.2 0.37 0.02 0.74 2.08 
300 m Regression Tal 0.02 —O.15 0.84 0.70 
Similancy $2.1 0.13 0.65 0.90 0.34 
id Gaussian 97.6 0.35 025 0.49 5.7 


Table 6. Evaluation of predictions of source emission rate, Q,(gs~‘), for ail Prairie Grass runs, using 
observed cross-wind integrated concentrations, C?, on monitonng arcs at distances of 50, 200 and 800m. 
Three different models are used to calcuiate Q, = C2,'(C7/Q), 


Fraction —S 

Distance to C, g.-9, Correlation within 20-2 
observation ” Model Os(gs~') Q, r factor of 2 0,0. 

Observed 68.4 

50m Regression 78.4 0.14 0.94 0.91 9.23 
Similanty 61.7 —0.10 0.98 1.00 . 0.02 
Gaussian 106.0 0.55 0.74 0,88 0.55 
200 m Regression 81.2 0.19 0.80 0.94 0.21 
Similarity 71.6 0.05 0.89 1.00 0.04 
Gaussian . 88.1 9.29 0.46 0.83 0.43 
800m Regression . 66.1 —0.03 0.07 0.74 0.46 
Similariry 69.2 0.01 0.60 0.88 0.26 
Gaussian 68.4 0.00 —0.15 0.65 2.98 


in this case is 16.2. Since @, (day)—@, (night)=53 gs~', 
it is conciuded that a significant difference in Q values 
could be discerned by this model if AQ, =Q,(day)- 
Q,(night) drops as low as 53 x (2.04/16.2)=6.7 gs7'; 
that is, if AQ, is about 15% of Q,, Therefore, in these 
best of research-grade experiments, where there are 
about 20 daytime and 20 night-time tests, and where a 
dispersion model is fit to these same data, a day-night 
difference in @ of less than 15% of the mean would not 
be estimated to be significant by this procedure. 


Predictions of source emission rates for ail cases 


Putting aside the questicn of differences in daytime 
and night-time averages in @, it is possible to use the 
complete database to estimate the overail ability of the 
models to estimate the 39 individuaj source emission 
cates. The reiative bias, the correlation, the fraction of 
predictions within a factor of two of observations, and 
the normalized mean square errof for the three models 
are listed in Table § and 6, for maximum point 
concentrations and cross-wind integrated concentra- 
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ons, respectively. These performance measures are 
calculated using data from the 50 m, 200 m and 800 m 
arcs. 

In most cases, the similarity model shows the most 
accuracy and the Gaussian mcdei shows the least 
accuracy. Also, the use of observed cross-wind integ-. 
tated concentrations, C?, leads to better results than 
the use of observed point maximum concentrauons. 
The deterioration of model performance as distance 
increases can also be seen. Focusing on Table 6, for 
cross-wind integrated concentrations, it is seen that 
the most accurate model, the similarity model. pro- 
duces a mean relative bias with a magnitude of 0.10 or 
less on all three distance arcs. The correlation drops 
from 0.98 to 0.60 as distance increases from [0 to 
300 m, and the fraction within a factor of :wo drops 
from 1.00 to 0.88 over the same distance. The normal- 
ized mean-square-error increases from (6.02 to 0.26 
over those distances, implying that the root-mean- 
square-error (rmse} increases from about 15% to 50% 
of the mean. In contrast, the Gaussian model yields a 
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much larger rmse that increases from about 70°% to 
300% of the mean. 

The factor of 5 or 6 difference in cmse between 
the similarity and Gaussian modeis impiies that the 
similarity model can discern differences in source 
emission rates that are a factor of 5 or 6 less than 
the minimum differences in source emission rates disc- 
erned by the Gaussian model. By comparing these 
numbers between Tables 5 and 6. it is seen that the 
relative rmses using C, observations are about 50% 
to 100% larger than the relauve rmses using C? 
observations. Consequently the use of C3 data permits 
one to discern differences in source emission rate that 
are 50% to 100% smaller than the minimum differ- 
éfces in source emission rates discemmed by C, data 


IMPLICATIONS FOR FUTURE RESEARCH 


The uncertainties in estimating source emission 
fates were first investigated using the Prairie Grass 
database, since the source conditions were simplified 
{a single continuous non-buoyant point source near 
the ground), the source emission rate was closely 
monitored, the site was fat and uniform. comprehens- 
ive meteorological data were taken, and observations 
of ground-level concentrations were taken at many 
points along monitomng arcs at ive downwind dis- 
tances. This experiment represents the optimum data- 
base of its type, and has been used in numerous 
research programs on atmospheric turbulence and 
dispersion. For this reason, the uncertainties in any 
analysis procedure should be minimized at this site. 
Conversely, if these procedures were to be applied to 
other experiments at other sites, where conditions are 
not so steady, smooth, or carefully-observed, the 
uncertainties would be expected to be larger. 

In the future, we will test these procedures using less 
ideal ooservations from other field experiments. These 
will include U.S. Army field tests of fog oil generators, 
which are used for smoke obscuration purposes. The 
plumes from fog oil generators are more complicated 
than those in the Praine Grass expermments, since the 
fag oil plumes are characterized by significant 
Momentum and buoyancy fluxes. Furthermore, the 
supporting meteorological data are not as complete as 
at Prairie Grass. Future tests will also include experi- 
mentaj data obtained for plumes from tall power piant 
stacks, where plume nse and mixing depth are com- 
plicating factors. From each set of experiments, the 
relative uncertainty of the procedures for estimating 
source emission rate mi] be assessed. It is expected 
that the uncertainty will increase as source conditions 
become more compiex or as input data are less 
complete. 

This study has demonstrated that models that do 
weil (in a least squares sense) in estimating point 
concentrations do not necessanly do weil when they 
are ‘turned around’ ‘to estimate source emission rate. 
The reason for this difference is the strong vanation 
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with distance of the plume centerline concentration. 
Consequently it is possible that a mode! which Aas 
zero mean dias in its concentration estimates. wiil 
have 2 50° mean bias in its source emission rate 
estimates. For the same reason. if a variable is Arst 
made non-dimensional and the model is ‘tuned’ with 
the data (e.g. the similamty mode! described above}, 
the predictions of the non-dimensional variable (e.g. 
Cu,x*/Q) may have zero bias, while the predictions of 
the concentration, C, may have significant bias. For 
optimum results, any tuning or regression analysis 
should be done with the variable that is of ultimace 
interest. 
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DISPLAY OF RELATIONS AMONG DATA USING BOX PLOTS 


One way of displaying large numbers of data is through the use of box 
plots, where the cumulative distribution function (cdf) of the dependent 
variables within a group of data is represented by a set of significant 
percentile values. For example, the 2th, 16th, 50th, 84th, and 98th 
percentiles are used in our analyses. These five significant points in the 
cdf are plotted by the SIGPLOT program using a “box" pattern as seen in the 
examples below. Variations of one type of data with another can be seen by 
breaking up the first type of data into groups defined by ranges of the second 


type of data. 


The SIGPLOT plotting package (see Appendix E) is used to generate the box 
plots. The ANADISTR program, described below, is used to generate the special 
input file required by SIGPLOT from a file containing multiple columns of 
data, representing concurrent values of variables such as observations of 
concentrations, wind speed, or stability. In the ANADISTR program, the user 
defines certain ranges. of the primary variables in the input file to be used 
for grouping the dependent variables and plotting them by means of box plots. 
The ANADISTR program requires one input file and generates one output file. 
The output file then serves as the input file to the SIGPLOT plotting package. 
There are no. default names associated with these files, and the user is 
prompted for the file names during the execution of the program. The ANADISTR 


program is written in FORTRAN 77. 


The input data file of the ANADISTR program could contain multiple 
columns of dependent variables (such as concurrent values of concentration 
observation) and other primary variables such as wind speed and stability. 
The ranges of the primary variables are also defined in the input file, to be 
used for grouping the data prior to generating the box plots. Table D-1 
describes the format of the input file. Figure D-1 shows an example of an 
imput file. Note that the ANADISTR program makes no corrections or 


substitutions for missing data; it is the responsibility of the user to 


provide valid data at each position. 


The output file of the ANADISTR program contains distributions (the 2th, 
16th, SOth, 84th, and 98th percentiles of the cdf) of the first variable as a 
function of the second variable. The information stored in this output file 
can then be plotted using the SIGPLOT plotting package (see Figure D-2 for an 


example). 


During the execution of the ANADISTR program, the following questions 
will be asked: 


. Name of the input file: 


The user must: specify the name of the input data file here. There is no 


default answer. 


. Name of the output file: 


The’ user must specify the name of the output file. There is no default 


answer. 


° An input file typically contains several columns of dependent and 
independent data. The ANADISTR program handles one such column or the 
ratio of any two columns of dependent data, specified by the user, at a 
time. This is accomplished by asking the user to select any two columns 
between 0 and MM, (see Table D-2). The distribution of the ratio of the 
numbers in these columns will be analyzed. Note that column "0" is 
simply all 1’s, and is not part of the input data file. Therefore, if 
the user wants to investigate the distribution of the dependent data in 
column 2, then two integers, 2 and O should be entered. If the user 
wants to investigate the distribution of the ratio of the dependent data 


in column 2 to the data in column 1, then 2 and 1 should be entered. 
° Implement a lower threshold for the dependent variable? (y/n): 


The user has the option of specifying a lower threshold for the one 
dependent variable or the ratio of the two dependent variables chosen 
above. This is sometimes necessary if the logarithmic scule is to be 
used and there are zero or minute values in the data whose distribution 
is to be analyzed. The default (i.e., hitting the RETURN key) answer is 


y". 


TABLE D-1. FORMAT OF THE MANDATORY INPUT DATA FILE OF THE ANADISTR PROGRAM. 
THE FOLLOWING KEY LETTERS ARE USED IN THE FORMAT COLUMN - FF: FREE 
FORMAT, C: CHARACTER, I: INTEGER, AND R: REAL. 


LINE NO. FORMAT DESCRIPTION 


1 FF/I . There are four integer constants in this line, 
representing the total number of observations (NN, < 
501), the total number of dependent variable (MM, < 
16), the total number of blocks (KK), and the total 
number of primary variables (NVAR, < 11). Note that 
KK is not used by ANADISTR since the blocking of data 
is performed internally according to the defined 
ranges of the primary variables. The limits on NN, 
MM, and NVAR are assigned in the program using the 
PARAMETER statements, and can be easily changed. 


2 FF/I There are KK integer constants in this line, 
representing the number of pieces in each block. The 
sum of all these integers should equal NN. Note that 
the information in this line is currently not used by 
the ANADISTR program. 


3 FF/I There are MM character constants in this line; each 
one can be at most eight characters long, containing 
the name of each of the dependent variables. All 


character constants must be enclosed in apostrophes. 


4 FF/T There should be KK character constants in this line. 
Each one can be at most 20 characters long, . 
containing the name of each of the blocks. All 
character constants must be enclosed in apostrophes. 
Note that the information in this line is not used by 
the ANADISTR program. 


LINE NO. FORMAT 
Next NN lines: 


FF/R 


Next NVAR lines: 
“EF/ 
I,C,R 


DESCRIPTION 


There are MM+NVAR real numbers in each line, with the 
first MM numbers representing the dependent 
variables, and the following NVAR numbers 


representing of the primary variables. 


Each line describes the way each of the NVAR primary 
variables is to be blocked. The first parameter is 
an integer (IXR, < 21), representing the number of 
ranges for the primary variable. The second 


parameter is a character constant, at most 40 


' characters long, enclosed in apostrophes, 


representing the name of the primary variable. The 
next [XR+1 real numbers, in numerical ascending 
order, define the boundaries of the ranges. For 
example, the following line: 

4 'u(m/s)’ 0. 2 5. 10. 20. 

means that wind speeds should be divided into four 
groups where the distribution of the dependent 
variables within each group is to be calculated. The 
first group is for those data when wind speeds are 
between 0. and 2. m/s, the second group is for wind 
speeds between 2. and 5. m/s, etc. 

The Limits on [XR are assigned in the program using 
the PARAMETER statement, and can be easily changed. 
Note that the sequence of the NVAR lines must be 
consistent with that of the last NVAR columns 
described in the previous section. As an example, if 
the MM+ith column in the previous section contains 
information for wind speeds, then the first line in 
this section should also contain grouping information 


for wind speeds. 


79 42 4 

39 40 
‘DEPVAR-A’ ’OEPVAR-3’ ‘OEPVAR-C’ ’ DEPVAR-D’ 
*SUBSET1’ ‘SUBSET2’ 


616.9 708.7 594.7 $16.5 11 3.9 soc. 2 
604.1 689.2 585.3 496.7 12 3.4 1000, 2 
868.0 674.8 580.3 516.8 13 3.5 1100. 2 
498.6 668.8 6§2.1 $48.3 a 3.8 1200. 2 
393.1 360.2 704.7 $81.3 1s 4.7 1300, 2 
409.0 740.9 570.1 621.4 16 $.2 1000. 3 
640.2 249.6 $10.1 553.5 17 S.4 1100. 3 
265.3 259.6 463.4 446.0 18 4.9 1100. 4 
192.7 91.6 131.0 485.0 19 4.2 1100. § 
1149.1 1217.5 1116.1 520.6 10 2.6 1600, 2 
972.8 1278.8 1178.1 536.9 1 3.2 1900. 2 
1137.5 1225.7 1081.7 617.4 12 3.8 1600. 2 
669.5 1052.8 905.1 637.3 13 4.5 1600. 2 
595.5 862.0 862.0 664.1 14 $.0 1900. 2 
741.2 $89.5 767.0 665.3 15 §.1 1860. 2 
612.6 602.4 728.2 672.4 16 $.0 1500. 3 
312.0 398.9 657.5 659.5 1? $.2 1500. 3 
400.2 340.2 412.3 586.0 1s $.1 1500. 4 
264.7 612.1 774.2 705.9 16 $.7 1400, 3 
290.0 428.4 787.3 708.8 17 S.1 1800. 3 
439.5 395.9 $12.3 602.4 138 S.1 2000. 4 
444.0 226.9 441.4 681.2 19 4.4 2000, 5 
178.1 216.6 486.1 925.4 20 4.6 2000. § 
202.3 126.2 288.6 522.9 21 4.9 2000. 6 
128.8 16.5 0.5 834.9 22 4.6 Qa. 6 
200.2 301.9 208.9 728.0 23 3.4 Q. 6 
358.3 481.3 3594.9 742.4 24 $.4 QO. 6 
611.1 1010.2 987.1 679.9 14 4.4 1500, 2 
499. 782.5 921.6 725.7 1S $.0 1500, 2 
$37.8 724.0 826.8 675.9 16 4.7 1800. 3 
220.0 $23.3 908.2 640.3 17 3.9 1800, 3 
479.2 337.5 788.6 $44.7 18 4.2 2000, 4 
133.2 198.3 383.1 738.5 19 3.1 1800, § 
98.2 167.3 213.5 1064.3 20 3.2 1800. 6 
92.5 104.6 142.2 741.2 21 3.1 1200. 6 
22.0 127.4 176.3 805.2 22 3.3. 1200. 6 
353.0 307.8 167.1 576.9 20 3.8 2000, 5 
338.0 280.9 188.4 225.3 22 2.3 2000. 4 
233.3 3585.39 234.9 719.2 22 2.4 2000 3 
198.3 32.7 184.9 745.2 23 3.6 2000. 6 
307.2 3.9 126.3 564.3 24 3.5 2000. 6 
313.7 2.9 9.3 467.2 L q22 3. 4 
265.2 2.9 2.9 703.9 2 3.6 2. 6 
298.6 329.9 484.6 695.3 4 $.2 a. 6 
$27.7 308.0 295.9 775.0 $ 4.7 Qa. 6 
454.1 301.0 1.0 995.6 6 2.9 Q. 6 
240.39 417.5 361.1 933.8 7 3.4 QO. 6 
$90.8 $79.3 144.2 666.5 a 3.1 1500. 5 
638.3 756.6 608.9 400.1 9 3.4 1800, 4 
949.8 1004.2 8065.4 $28.9 10 3.4 1800. 3 
886.8 8385.6 706.2 $17.4 11 3.0 1300. 2 
635.5 761.0 670.9 $96.6 12 4.5 1200. 2 
359.3 412.6 232.5 937.6 1 2.3° 1200. 6 
484.7 360.7 226.8 973.0 2.5 1200. 6 
$29.7 332.0 202.5 980.9 3 2.4 1200. § 
$85.8 291.4 i86.1 1100.1 4 2.1 1200. § 
367.7 368.0 260.2 1005.5 5 2.1 1200. 6 
324.7 270.9 72.7 1088.6 6 2.0 1200. 6 
489.0 274.6 208.5 942.2 7 2.6 1200. 6 
$70.8 337.1 218.9 646.5 8 2.8 1200. 3§$ 
419.7 294.4 206.1 344.9 3 4.3 1200. 4 
532.8 414.2 197.9 477.9 9 4.8 1800. 3 
425.2 365.7 198.7 469.5 10 7.1 #1700. 4 
467.5 411.5 228.5 455.3 ll 7.3 2000. 4 
362.2 306.4 147.6 405.2 12 $.1 2000. 4 
429.2 287.4 139.2 480.6 13 5.4 2000. 4 
446.0 338.1 169.5 461.2 14 5.7 2000. 4 
292.9 253.8 148.6 460,7 15 $.8 2400. 4 
$30.3 322.5 257.2 460.5 16 7.3. 2700, 4 
364.9 326.7 251.1 $10.6 17 7,8 #43000, 4 
111.4 196.4 248.5 0.0 23 1.9 280. 4 
89.8 146.5 254.9 0.0 24 1.9 250, 4 
82.5 248.0 160.9 0.0 1 2.9 250, 4 
296.5 233.2 193.2 0.90 2 2.8 250. 4 
218.4 299.7 165.0 0.0 3 2.9 250. 4 
484.5 274.2 184.0 0.0 4 3.5 250. 4 
384.7 324.6 163.2 9.0 s 365 280. 4 
253.2 488.3 175.6 0.0 6 358 280. 4 
289.5 304.1 193.1 a.9 Agt 280. 4 
6 ‘nour of day’ -0.01, 4., 8., 12., 16., 20., 24.91 
10 ‘yu (m/si’ 0.5, 1.3, 2.5, 3.5, 4.5, 3.35, 6.5. 7.5, 9.5, 9.8, 10.5 
6 ’h (m)’ -0.01, 200., 600., 1000., 1500., 2000., 3000.: 
3 ‘pq class’ 0.5 3.5 4.5 6.5 


Figure D-1. An example of the input data file for the ANADISTR program. 
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DEMONSTRATION OF THE RESULTS GENERATED BY THE ANADISTR PROGRAM. 
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An example of the results generated by the ANADISTR program and 
plotted using SIGPLOT. Significant points on each box represent 
the 2nd, 16th, S5SOth, 84th, and 98th percentiles. The dashed 


lines represent the factor of two lines. 
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The next question will be asked only if the user answers "“y" to the 


above question. 


Enter the lower threshold (e.g. 0.01): 


There is no default answer, but 0.01 has proven to be a good choice for 


the ratio of two dependent variables in our tests. 


Implement an upper threshold for the dependent variable? (y/n): 
This is sometimes necessary if the logarithmic scale is to be used and 
there are very large values in the data whose distribution is to be 


analyzed. The default (i.e., hitting the RETURN key) answer is "y". 


Thé following question will be asked only if the user answers "“y" to the 


above question. 


Enter the upper threshold of the ratio (e.g. 100.): 


, There is no default answer, but 100 has proven to be a good choice for 
the ratio of two dependent variables in our tests. 


Intentionally Blank 
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APPENDIX E 


A. USER'S GUIDE FOR THE SIGPLOT PLOTTING PACKAGE 


The SIGPLOT plotting package developed at Sigma Research Corporation is a 
versatile tool for producing different kinds of two-dimensional plots, such as 
scatter plots, graphs, box plots (sometimes called residual or whisker plots), 
or error bar plots. The user can specify many parameters including the number 
of frames per page, the aspect ratio of the frame, and the mapping of the 
coordinates. The graphics library routines used by SIGPLOT, together with the 
screen and printer drivers (described later) were originally developed by Dr. 


¢ Arlindo daSilva of the University of Wisconsin at Milwaukee. 


¢ SIGPLOT requires two input files: 1) the template file that contains the 
control parameters which influence the appearance of the plots, and 2) the input 
data file that contains the data to be plotted. Tables E-1 and E-2 describe the 
formats of the template file aid the input data file, respectively. Examples of 
the template file are shown in Figures E-! and E-2 Examples of the input data 


file are shown in Figures E-3 and =-4. 


SIGPLOT creates a Tektronix picture file that can be viewed directly on 
any kind of the PC -raphics environments (e.g., Hercules, CGA, EGA, and VGA) 
using the screen driver, TEXPC. Hard copy output can also be generated from 
the Tektronix picture file with a printer driver. There are three printer 
drives, TEXEPS, TEKELQ,’ and PS, that are currently available. The first two 
drivers are used to drive an EPSON-compatible dot matrix printer, with TEXEPS 
for low resolution and TEKELQ for high resolution. The PS program is used to 
drive a PostScript printer, such as Apple LaserWriter, NEC LC-890, or TI 
MicroLaser PS35. It is recommended that the user have access to a PostScript 


printer to obtain the best results in the shortest time. 


* SIGPLOT requires about 200KB of memory. The other screen and printer 
drivers require less than 100KB of memory, except for TEXELQ, where 450KB of 
memory is required due to the high resolution and the use of the bitmap 
approach in the driver program. The SIGPLOT plotting package and the graphics 
library routines were written in FORTRAN. The screen and printer drivers were 


written in C. 


TABLE E-1. THE FORMAT OF THE TEMPLATE FILE OF SIGPLOT. THE FOLLOWING KEY 
LETTERS ARE USED IN THE FORMAT COLUMN - FF: FREE FORMAT, 
C: CHARACTER, I: INTEGER, AND R: REAL. 


The global control parameters are specified in the first section of the 
template file, lines 1 through 16. 


LINE NO. FORMAT . DESCRIPTION 

1-3 . Reserved for comments 

4 FF/C Name of the input data file, currently not used 
S FF/C Name of the output Tektronix picture file, 


currently not used 


6 FF/I Flag for the frame aspect ratio, 1-5, 
Ls: xy = lil 
2: xy = 1:2 
3: x:y = 2:1 
4: xy = 1:3 
8: xy 23:1 
7 FF/T : Number of frames per page, 1-4 
38-9 A80° Title for the page (no title will be drawn if 


“O" appears as the first character of the line) 


10 FF/C Flag (PAXIS) for the axis along which the first 
, column, representing the independent variable, 
of the data in the input data file (see Table 
E-2) will be plotted (x or y). PAXIS must = x 
if JPATIN (described below) = 4, and PAXIS must 
= y if IPATIN = 6 


11 FF/I Flag for mapping, 1-4, 
1: linear in x, linear in y 
2 linear in x, logarithmic in y 
3: logarithmic in x, linear in y 
4 logarithmic in x, logarithmic in y 


E-2 


TABLE E-1. 


LINE NO. 


12 


13 


i4 


is 


16 


THE FORMAT OF THE TEMPLATE FILE OF SIGPLOT. THE FOLLOWING KEY 
LETTERS ARE USED IN THE FORMAT COLUMN - FF: FREE FORMAT, 


C: CHARACTER, 


FORMAT 


FF/T 


FF/I 


FF/C 


SA1 


FF/I 


I: 


INTEGER, AND R: REAL. Continued. 


DESCRIPTION 


Flag (IPATIN) for plot pattern, 1-6, 


1: scatter plot 
2: line graph 
3: scatter plot except line graph for the 


last variable 
box plot 
error bar plot 


same as S but with extra labelling 


Flag for background, 0 or 2, 
0: no background 
2: gridded background 


Flag for system time, y or n, 
if y: system time will be printed out on the 


upper right corner of each page 
Five point patterns for the ‘scatter plot , 


Flag (IEXTRA) for the plotting of extra lines, 


I: , x=0 will be plotted 

2: y=0 will be plotted 

3: x=0 and y=0 will be plotted 

4: x=1l will be plotted 

3: y=1 will be plotted 

6: x=1 and y=1 will be plotted 

7: diagonal line will be plotted 

8: y=0.5 and y=2 (factor of two) will be 
plotted 

9: x=-0.667, 0, and 0.667, and 


ye4x2/ (4x7). (see text) will be plotted, 


else: no extra lines will be plotted. 


Note that IEXTRA = 9 is effective only 
if IPATIN = 6 
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TABLE E-1. THE FORMAT OF THE TEMPLATE FILE OF SIGPLOT. THE FOLLOWING KEY 
LETTERS ARE USED IN THE FORMAT COLUMN - FF: FREE FORMAT, 
C: CHARACTER, I: INTEGER, AND R: REAL. Continued. 


The next section of the template file (lines 17 through 29) contains the 
parameters that are applicable to a frame. This section can be repeated if 
there are multiple frames to be plotted in a print job. However, the use: can 
prepare just one such section if the same information is to be used repeatedly 


by all frames. 


LINE NO. FORMAT DESCRIPTION 

17-19 Reserved for comments 

20 FF/R Constants, a and b, for the linear 
transformation of the independent variable, 
where 
*new ~ birt mies 


a=1 and b=0 means no transformation is needed 


21 FF/R Constants, a and b, for the linear 
transformation of the first dependent variable, 
where 


Y1,new . 3") old eras 


a=1 and b=0 means no transformation is needed 


22 FF/R. Same as above, but for the second dependent 


variable 

23 FF/R Same as above, but for the third dependent 
variable 

24 FF/R Same as above, but for the fourth dependent 
variable 
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TABLE E-1. THE FORMAT OF THE TEMPLATE FILE OF SIGPLOT. THE FOLLOWING KEY 
LETTERS ARE USED IN THE FORMAT COLUMN - FF: FREE FORMAT, 
C: CHARACTER, I: INTEGER, AND R: REAL. Concluded. 


LINE NO. FORMAT DESCRIPTION 


25 FF/R Same as above, but for the fifth dependent 
variable. Note that lines 22 through 25 cannot 
be omitted even if only one group of data were 
to be plotted 


26 FF/R xmin, xmax, and dx of the x-axis 

27 FF/R ymin, ymax, and dy of the y-axis 

28 FF/C Format specifier for the numerical labels of the 
x-axis. If "!" appears as the first character of 


the line, the appropriate format will be 
determined internally by the program; otherwise, 
the user should supply a simple FORTRAN [-, F-, 
or E~format specifier, enclosed in parentheses, 
e.g., (15), (F6.3), and (E8.1) are accepted, but 
(315), (15,f6.3), (1P,£8.1), and (G9.1) are not 


accepted 


29 FF/C Format specifier for the numerical labels of the 


y-axis 


TABLE E-2. THE FORMAT OF THE INPUT DATA FILE OF SIGPLOT. ‘THE FOLLOWING KEY 
LETTERS ARE USED IN THE FORMAT COLUMN - FF: FREE FORMAT, 


C: CHARACTER, 
LINE NO. FORMAT 
1 A40 
2 A40 
3 A40 
4 FF/I 


Next NPTS Lines: 


For [PATIN = 1, 


2, and 3, 


FF/R 


I: INTEGER, AND R: REAL. 
DESCRIPTION 


Title for the frame (no title will be drawn if “O" 


appears as the first character of the line) 


Label for the x-axis (no label will be drawn if "0" 


appears as the first character of the line) 


Label for the y-axis (no label will be drawn if “O" 
appears as the first character of the line) 

Two integers specifying the number of points (NPTS) 

and the number of groups of data (MANY) to be plotted. 
MANY cannot be > 5 for IPATIN = 1, 2, 3, and 5, and 
MANY must be = 1 for IPATIN = 4 and 6. NPTS cannot be 

> 700 for IPATIN = 1, 2, and 3. NPTS cannot be > SO for 
IPATIN = 4, 5, and 6 (see text). 


There are 1+MANY real numbers in each line. The first 
number represents the independent variable, which can be 
plotted either along the x- or the y-axis depending the 
value of PAXIS (see Table E~1). The next MANY numbers 
represent the dependent variables. For example, if 

three curves (MANY=3), £, (x), 5 (x), and f(x) were to be 
plotted, then each line here should contain four real 

Li foi and 4 i where i=l, NPTS. If 
PAXIS = "x", the x will be plotted along the abscissa, 
and f. f,. and f, 


vice versa PAXIS = 


numbers, Kyo £ 


will be plotted along the ordinate; 


y 


TABLE E-2. THE FORMAT OF THE INPUT DATA FILE OF SIGPLOT. THE FOLLOWING KEY 
LETTERS ARE USED IN THE FORMAT COLUMN - FF: FREE FORMAT, 
C: CHARACTER, I: INTEGER, AND R: REAL. Continued. 


LINE NO. FORMAT DESCRIPTION 
For IPATIN = 4, 


FF/R, I There are six real numbers and one integer in each line. 
The first real number represents the independent 
variable. The next five real numbers rerresent the 
values of the dependent variable at the 2th, 16th, lvth, 
84th, and 98th percentiles, respectively. Note that the 
value of the independent variable listed here frequently 
represents a range of the independent variable; for 
example, a wind speed of 7 m/s actually represents wind 
speeds in the range of 5 to 8 m/s. The integer 
represents the number of data points based on which the 
distribution of the dependent variable is derived. No 
box will be plotted if the number of data points is less 
than five since not enough information is available to 
define a distribution. 


For IPATIN = 5, 


FF/R There are 1+3#MANY real numbers in each line. The first 
humber represents the independent variable. The 
remaining numbers for the dependent variables are in 
MANY groups of three numbers. The three numbers, which 
must be in order, represent the distribution of a 
dependent variable. This distribution can be 1) u-c, x, 
and ut+o, where uw is the mean, and o is the standard 
deviation, or 2) lower c.1l., nominal value, and upper 


c.1l., where c.l. is the confidence limit. 
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‘TABLE E-2. THE FORMAT OF THE INPUT DATA FILE OF SIGPLOT. THE FOLLOWING KEY 
LETTERS ARE USED IN THE FORMAT COLUMN - FF: FREE FORMAT, 
C: CHARACTER, I: INTEGER, AND R: REAL. Concluded 


LINE NO. FORMAT DESCRIPTION 


For IPATIN = 6, 


FF/R,C There are four real numbers and one character constant 
(no more than 17 characters long) in each line. The 
definition of the first four real numbers is identical 
to that when IPATTN = 5, except now MANY must = 1. The 


character constant, enclosed in apostrophes, is used to 
label each data point. 


The above 4+NPTS lines provide enough information to plot a frame. Additional 
data, similar in structure, can be appended here if the plotting of more than 
one frames in a print job is desired. 


ee es he er ee oe re eee ey ey ee ee Oe a 


! 
! Main switches for plotting. 


t ce} om 2 Q enn nn =- we Q ee wee nn enw e enw een n ence wooo sens 
urrs.l Name of input data file. 

tekl.pic Name of oucput tektronix file. 

1 Aspect ratio (integer, 1 - 5). 

Bt Number of plots per page (integer, 1 - 4). 

demo of ipattn=2 

) 

x Which axis serves as independent variable (x or y). 

1 Flag indicating log or linear mapping (1 = 4). 

2 Pattern. 

0 Background specification. é 

y Print out system time on the upper right hand corner (y or no). 
«toes Patterns of scatter plots (Sal) 

0 Extra line,1:x=0,2:y=0,3:x, yw0,4:xel,S:y=l, 6:x,yel, 7:diag, 8:y=fac. 2.,9:f>-nmse, else:nothing. 
§ wm 288 Oe BRST oeewseooenesonoe ee BOBS SSOP esseseenasasoansoenacun 
H Parameters for plot 1. 

{ mee: _ PE SD CS PN ED OP SO DES OE AY OP OS OO IPD OP OD OP OE A EE DE OSMOSIS BEMsEBewwawe =e 
1. 0. ascale, bacale for the independent variable axis. 

1. 0. ascale, bscale for curve 1. 

1. 0. ascale, bacale for curve 2. 

1. 0. ascale, bscale for curve 3. 

1. 0. ascale, becale for curve 4. Z 
1. 0. ascale, bscale for curve 5. 


-6.28319 6.28319 3.141595. xmin, «max, and dx for the x axis. 
“1.2 1.2 0.3 ymin, ymax, and dy for the y axis. 

(£8.2) format. for x label 

(£4.1) format for y label 


Figure E-1. ‘An example of the template file of SIGPLOT. Refer to Figure E-6 for 


the results. 


lecwesccncneuen sada acansccewsekawaasatwacubuawensseveuscesecenwcuucescous 


! Main switches for plotting. 


! Qo. x1} Yemen wane an meee nee nen eee me eons 
urrs.1 Name of input data file. 

texl.pic Name of output tektronix file. . 
1 Aspect ratio (integer, 1 - 5). 

1 Number of plots per page (integers, 1 - 4). 

demo of tpattne§ 

9 

x Which axis serves as independent variable (x or y). 

3 Flag indicating log or linear mapping (1 =- 4). 

A) Pattern. 

fe) Background specification. (0 or 2) 

Fy Print out system cime on che upper tight hand corner (y or a). 
70.aS Patterns of scatter plots (Sal) 

QO Extra line, i:x#0,2:y00,3:x,y=0,4:xel,S:yel, 6:x, yal, 7:diag,8:ysfac. 2.,3:So-nmse, else:nothiag. 
Y even an an cnepereenenananaaran anes aves: we we eeaenenan. pen an Gresesapenan: aaweawasma 
{ Parameters for plot 1. 

! a a as a <7 ce an a ED A 
1. 0. ascale, bacale for the independent variable axis. 

1. 0. ascale, bscale for curve 1. . 

1. 0. aacale, bscale for curve 2. 

Le ue ascale, bacale for curve 3. 

2. 0. ascale, bscale for curve 4. 

1. 0. ascale, dscale for curve $. 


200. 20000. 10. xmin, «max, and dx for the x axis. 
“1.5 1.5 0.5 ymin, ymax, and dy for the y axis. 
(45) 

(£4.23) 


Figure E-2. An example of the template file of SIGPLOT. Refer to Figure E-9 for 


the results. 
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50 5 


-6.03186 
-5.78053 
+5 .52920 
-5.27788 
~§.02655 
74.77822 
=4.52389 
=4.27257 
4.92224 
~3.76991 
-3.51858 
@3.26726 
-3.01593 
~2.76460 
72.51327 
=2.26195 
-2.01062 
-1.78929 
-1.50796 
“1.25664 
-1.00531 
-0.753983 
70.502655 
“0.251328 
3.9090C9 
0.251328 
0.502655 
0.753982 
1.00531 
1.25664 
1.50796 
1.75929 
2.01062 
2.26195 
2.51327 
2.76460 
3.91593 
3.26726 
3.51858 
3.76991 
4.92124 
4.27257 
4.52389 
4.77522 
$.02655 
$.27787 
$.$2920 
$.78083 
6.03186 
6.28319 


Figure E-3. 


0.248690 
0.481754 
0.684547 
0.844328 
0.981087 
0.998027 
0.982287 
0.904827 
0.770823 
0.587785 
0.368125 
0.228333 
-0.125333 
-0.368125 
-0.587785 
-0.770513 
-0.904827 
0.982287 
-0.998027 
-0.981057 
-0.844328 
-0.684547 
-0.481754 
~9.248690 
9.300000 
2.248690 
0.481753 
0.684547 
0.844328 
0.981056 
0.998027 
0.982287 
0.904827 
0.770813 
0.587785 
0.368124 
0.325334 
-0.125333 
-0.368124 
~0.587785 
-0.770513 
0.904827 
-0.982287 
-0.998027 
-0.951057 
~0.944328 
-0.684548 
~0.481754 
0.248690 
~0.301992E-06 


An example of the input data file of SIGPLOT. 


for the results. 


0.368124 -0.844328 
-0.125333 -0. 684547 
0.125333 0.481753 
0.368125 0.248690 
0.587786 0.397359E-06 
0.770513 0.248690 
0.904827 0.481754 
0.982287 0.684547 
0.998027 0.344328 
0.951086 0.951087 
0.844328 0.998027 
0.684547 0.382287 
0.481753 0.904827 
0.248690 0.770513 
9.397391E-07 0.587785 
-0.248690 0.368124 
0.481754 0.128333 
0.684547 -0.128334 
-0.844328 0.368124 
-0.951057 -0.587785 
0.998027 -0.7708133 
-0.382287 -0.904827 
-0,904827 -0.982287 
9.770513 -3.398027 
-0.587785 -0.951056 
-0.368124 -0.844328 
-0.125333 -0.684547 
0.125333 0.481754 
0.368125 ~0.248690 
0.587785 -0.381470E-06 
0.770513 0.248690 
0.904827 0.481754 
0.982287 0.684547 
0.998027 0.844328 
0.951056 0.951057 
9.844328 0.998027 
0.684548 0.982287 
0.481754 0.904827 
0.248690 0.770514 
O.341731E-06 0.587785 
0.248690 0.368125 
-0.481754 0.125333 
-0.684547 0.128333 
-0.844327 0.368124 
0.951056 -0.587785 
0.998027 -0.770813 
-0.982287 0.904627 
0.904827 -0,982287 
-0.770513 -0.998627 
0.587785 -0.951087 
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0.398027 
0.382287 
79.904827 
70.770513 
0.587785 
-0.368125 
0.125333 
0.325333 
0.368125 
0.587785 
0.770813 
0.904827 
0.382287 
0.998027 
0.951056 
0.844328 
0.684547 
0.481753 
0.248690 
7-0.§56284E~-07 
-0.248690 
79.481753 
79.684547 
“0.344323 
-9.98i557 
70.998027 
0.382287 
0.904827 
-0.770513 
0.587796 
0.368125 
-9.128333 
0.125333 
0.368125 
0.587785 
0.770813 
9.904827 
0.982287 
0.998027 
0.951057 
9.844328 
9.684547 
0.481754 
0.248691 
0.723200E-06 
-9.249689 
9.481753 
0.684847 
9.844328 
0.951056 


0.770814 
“9.304827 
9.982287 
0.998027 
79.951056 
70.844328 
0.684547 
0.481753 
~0.248690 
0.2543C8E-06 
9.248690 
0.481754 
0.684547 
0.344328 
0.951957 
0.998027 
0.982287 
0.304827 
0.770513 
0.587788 
3.368124 


9.770533 
~0.904827 
0.982287 
-9.998027 
0.951057 
0.844328 
~0.684547 
70.481754 
~0.248690 
0..90735E-06 
9.248690 
9.481753 
2.684547 
2,844328 
9.951056 
2.998027 
0.982297 
9.904827 
0.770514 
0.587736 
0.368125 
0.425334 
~0.325333 
-0.368128 
-0.587785 


Refer to Figure £-6 


all periods 
n-a distance (m) 


var (dws) 
6 3 
312 
625 


/ median 


«5 -0.029 
20 -0.009 


1250.0 -0.035 


2500 
$000 
10000 


Figure E-4. 


0 0.155 
«8 -0.033 
«0 -0.417 


lemin var (ws) 


0.002 
0.005 
0.012 
0.015 
6.170 
0.061 


0.034 
0.019 
0.059 
6.188 
6.373 
0.539 


0.943 
0.048 
0.070 
0.100 
0.323 
0.369 


0.112 
0.085 
0.148 
0.325 
0.633 
0.876 


0.166 
0.209 
0.330 
0.359 
0.440 
0.406 


Q.275 
0.361 
0.480 
0.565 
0.752 
0.873 


An example of the input data file of SIGPLOT. 


for the results. 
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0.383 
0.513 
0.630 
0.771 
1.062 
1.339 


Refer to Figure E-9 


AS one can see from Table E-1, SIGPLOT is capable of creating the 


following kinds of plots: 


IPATIN = 1: scatter plot (e.g., Fig. E-5) 

IPATIN = 2: line graph (e.g., Fig. E-6) 

IPATIN = 3: scatter plot except line graph for the last variable 
(e.g., Fig. E~7) 

IPATIN = 4: box plot (e.g., Fig. E-8) 

IPATIN = 5S: error bar plot (e.g., Fig. E-9) 

IPATIN = 6: same as IPATIN = 5 but with extra labelling (e.g., Fig. 


E-10) 


The usage of each option is described below. 


For IPATIN = 1, groups of data are represented by different dot patterns 
that are defined in the template file (see Table E-i). At most, five grours 
ef data (MANY = 5) can be plotted, with a maximum of 700 points for each 


group. 


IPATIN = 2 is similar to IPATTN = 1 except that points are now connected. 
The following line patterns are used to represent different curves: solid, 
short-dashed, long-dashed, dot-dashed, and dotted. At most, five curves 
(MANY = 5) can be plotted, with a maximum of 700 points for each group. No 
user customization of the line patterns is allowed. I[t is important that the 
data points in the input file are sorted according to the independent 


variable. 


IPATIN = 3, a combination of IPATIN = 1 and 2, is useful when the user 
wants to see how well -a theoretical curve fits the observed data. Although 
the order of the data points does not matter for a scatter plot, in this case 
it is important that the data points in the input file are serted according to 
the independent variable. At most, five groups of data (MANY = 5) can be 


plotted, with a maximum of 700 points for each group. 


The IPATIN = 4 option is an alternative to the scatter plot when the 


number of data points is large. In preparing the inppt data file for SIGPLOT, 


E-13 


the user first defines certain ranges of the independent variable to be used 


for grouping the dependent variable. The distribution of the dependent 
variables within each group is then determined and represented by five 
Significant points in the cumulative distribution function (cdf). These five 
values could be the 2th, 16th, SOth, 84th, ana 98th percentiles of the cdf, or 
the mean and mean + one and two standard deviations. SIGPLOT then uses a box 
pattern to represent the distribution of the dependent variable within each 
grouping or range of the independert variable. Only one set of data (i.e., 
MANY = 1, even though five points are needed to define a box) is accepted for 


this option, with a maximum of 50 boxes. 


IPATTIN = S is similar to IPATIN = 4 except that three values (vs. five) 
are needed to define an error bar (vs. a box). These three values can be the 
Mean and mean + one standard deviation of a dependent variable, or the nominal 
value of a dependent variable and its 95% confidence limits. At most, five 
groups (MANY = 5) of data can be plotted, with a maximum of 50 error bars for 
each group. The following error bar patterns are used: filled square, empty 


square, filled triangle, empty triangle, and cross. 


IPATIN = 6 is similar to IPATIN = 5 except that the user can label each 
data point. Because of the additional information to be plotted, only one 
group of data (MANY = 1) is accepted, with a maximum of 50 error bars. This 
option is designed primarily to plot the FB (fractional bias), together with 


its confidence limits, against the NMSE (normalized mean square error), where 


FB = P (E-1) 


(E-2) 


{If [EXTRA = 9 (see Table E-1), SIGPLOT will plot the additional 
x = -0.667, 0, and 0.667 lines, representing the factor of two and zero FB 
lines, together with the y = 4x*/(4-x7) line, representing the “minimum" NMSE 


(due only to the mean bias) as a function of FB. 


. 
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The information concerning the usage of the driver programs, TEXPC, 
TEKELQ, TEKEPS, and PS, can be obtained by simply executing the programs 


without providing any arguments, and will not be repeated here. 


Finally, an example is given below of the procedures followed to use the 


graphics package. 


Step 1: 


Step 2: 


Step 3: 


Step 4: 


Step 5: 


The user prepares the template file (DEMO.INQ) and the input data 
file (DEMO.DAT) according to the formats described in Tables E-i and 
E-2. The user can create his own template file by editing the 
sample template file. The input data file is usually generated by 


some other programs. 


After the execution of SIGPLOT, a Tektronix picture file (DEMO. PIC) 
is generated. 


The user can view the results on screen by typing: 

TEXPC DEMO. PIC 

if a Hercules graphics card is installed, or 

TEKPC DEMO.PIC 16 

if an EGA (with a resolution of 640x350 pixels) graphics card is 
installed. 


A high resolution hard copy output can be generated on an 
EPSON-compatible dot matrix printer by typing: 
TEKELQ DEMO. PIC. 


Or if the user has access to a PostScript printer, a PostScript file 
(DEMO. PS) will be created by typing: 

PS DEMO. PIC, 

and this file can be printed out by typing: 

PRINT DEMO. PS 
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OEMO OF [PATTN=1 


URBAN SITE 1287, OBSERVED GH 


Us PREDICTED (H/S) 


68 8.2 84 8.6 8.8 1.8 1.2 
Us OBSERVED (M/S) 


URBAN SITE 187. PREDICTED GH 


Us PREDICTED (H/S) 
a 
o 


a8 8.2 8.4 #85 @&8 1.8 1.2 
Us OESERVED (M/S) 


Figure E-5. A sample scatter plot (IPATIN = 1) generated by SIGPLOT. 


E~16 


DEMO OF IPATTN=2 


- @3/25/91 


-6.28 ~3.14 8.22 3.14 6.28 

xX 
Figure E-6. A sample line graph (IPATIN = 2) generated by SIGPLOT. Refer to 
Figures E-1 and E-3 for the template and data files used for this 


figure. 
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SIGMA-V/Ux 


CEMO OF I[PATTN=3 


OBSERVED (DOTS) VS. THEORETICAL (LINE) 
180.2 


19.8 


3u2 


1.8 


6.3 


0.1 
Q.1 Q.3 1.8 3.2 18.0 31.6 188.0 


Wa /Use 


Figure E-7. A sample scatter plot and line graph (IPATIN = 3) generated by 


SIGPLOT. 
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DEMO OF IPATTN=4 


188.28 188.28 
18.28 13.28 
2 2 
ra) 
| 8 
~ 1.38 > 1.38 
: 8 3 
z = 
8.18 8.18 
8.21 8.81 
WINO SPEED (M/S) STABILITY INCEX 
188. 38 . 188.88 
18.38 18.28 
a) 7) 
8 8 
> 1.28 = 1.28 
~] ° 
z z 
a. 18 8.19 
8.81 8.31 


MIXING HEIGHT =198 (M) 


Figure E~8. A sample box plot (IPATTN = 4) generated by SIGPLOT. 
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DEMO OF IPATTN=S 


ALL PERIOOS 


VAR (DHS) / MEDIAN 1-MIN VAR(WS) /2 


200 2000 29008 
N-S DISTANCE (M) 


Figure E-9. A sample error bar plot (IPATIN = 5) generated by SIGPLOT. Refer to 
Figures E-2 and E-4 for the template and data files used for this 
figure. . 
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DEMO OF IPATTN=5 


THORNEY ISLAND. INSTANTANEOUS 
15.8 


izes 


12.8 


tes 


NMSE 


5.8 


29 


“2.8 -1.5 -1.8 -@.5 4.2 8.5 1.8 i) 2.8 
FB (WITH SS-PERCENT C.1.) 


Figure E-10. A sample error bar plot with labelling (IPATIN = 6) generated by 
SIGPLOT. : . 


“E-21 


Intentionally Blank 
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APPENDIX F~2 


LISTINGS OF THE DUGWAY DATA ARCHIVES - SMOKE/OBSCURANT DATASETS 
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